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Abstract: 

The  dynamics  of  waves  in  deep  water  and  in  shallow  water  are 
quite  different.  In  shallow  water,  we  have  shown  experimentally  the 
existence  of  a  family  of  finite-amplitude  waves  that  propagate 
practically  without  change  of  form  in  shallow  water  of  uniform 
depth.  The  surface  patterns  of  these  waves  are  genuinely  two- 
dimensional,  and  periodic.  The  waves  are  easy  to  generate 
experimentally,  and  they  seem  to  be  stable  to  perturbations.  They 
are  described  with  reasonable  accuracy  by  an  8 -parameter  family  of 
exact  solutions  of  the  Kadomtsev-Petviahsvili  equation. 

The  situation  is  quite  different  in  deep  water,  where  much  of 
our  knowledge  is  based  on  numerical  simulations.  An  approximate 
model  of  one-dimensional,  nearly  monochromatic  waves  in  deep 
water  is  the  nonlinear  Schrddinger  equation.  We  find  that  for  certain 
ranges  of  parameters  in  initial  data,  numerical  solutions  of  the 
equation  are  so  unstable  that  long-time  simulations  of  the  equation 
are  not  reproducible,  and  are  quite  unreliable.  The  equation  is 
integrable,  and  not  chaotic,  but  numerical  simulations  of  the  equation 
are  chaotic.  We  call  this  "numerically  induced  chaos". 
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Nonlinear  Ocean  Waves 
Final  report:  If 1192  -  9/30/94 

The  scientific  study  of  water  waves  dates  back  to  the  works  of 
Newton,  the  Bemoullis,  and  Euler,  but  the  observations  and 
speculations  of  ancient  mariners  long  preceded  these  works.  The 
mathematical  problem  of  invisicid  water  waves  was  first  posed 
properly  by  Stokes  (1845;  see  Stoker,  1957,  p.  16).  Stokes'  equations 
are  nonlinear,  but  Stokes  also  linearized  his  equations  and  found 
both  linear  and  weakly  nonlinear  approximate  solutions,  for  waves  of 
small  amplitude.  Nekrassov  (1921)  and  Levi-Civita  (1925)  showed 
that  Stokes'  approximate  solutions  represent  exact  solutions  of  the 
nonlinear  problem. 

Stokes'  small-amplitude  solutions  notwithstanding,  nonlinearity 
plays  an  essential  role  in  the  dynamics  of  almost  all  water  waves.  In 
shallow  water,  waves  have  wavelengths  longer  than  the  mean  depth 
of  the  fluid,  and  nonlinearity  leads  either  to  wavebreaking  (Airy, 
1845)  or  to  the  formation  of  solitary  waves  (Russell,  1838,  1844).  In 
deep  water,  waves  have  wavelengths  shorter  than  the  mean  depth  of 
the  fluid,  and  Benjamin  &  Feir  (1967)  showed  that  the  plane  waves 
studied  by  Stokes  and  others  are  subject  to  a  nonlinear  instability. 
Thus,  nonlinearity  plays  an  essential  role  in  the  formation  of  almost 
all  waves  that  are  observed  in  the  ocean,  no  matter  what  their 
wavelengths. 

During  the  nearly  three-year  duration  of  this  grant,  four 
principal  investigators  (Mark  Ablowitz,  James  Curry,  Joe  Hammack, 
Harvey  Segur)  and  their  colleagues  have  studied  aspects  of  nonlinear 
water  waves,  in  both  shallow  and  deep  water.  Phenomena  in  deep 
water  differ  significantly  from  those  in  shallow  water,  and  our 
studies  in  these  two  areas  are  distinct.  These  studies  combine 
mathematical  analysis,  laboratory  experiments,  numerical 
computation,  and  field  observations.  In  this  report,  we  summarize 
our  major  published  conclusions,  and  we  outline  briefly  the  current 
status  of  work  in  progress. 

A.  Shallow  water 
1.  Published  work 

The  focus  of  our  work  has  been  on  spatially  periodic  waves  of 
permanent  form  that  appear  to  be  stable  in  shallow  water.  A  full 
discussion  on  this  work,  which  will  appear  in  the  Journal  of  Fluid 
Mechanics,  can  be  found  in  Appendix  A  of  this  report.  Our  main 
conclusions  are  as  follows. 
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•  Laboratory  experiments  demonstrate  the  existence  of  a  family  of 
gravity-induced,  finite-amplitude  water  waves  that  propagate 
practically  without  change  of  form  in  shallow  water  of  uniform 
depth.  The  surface  patterns  of  these  waves  are  genuinely  two- 
dimensional,  and  periodic.  The  basic  template  of  a  wave  is 
hexagonal,  with  six  wave  crests  forming  the  boundary  of  a  cell  that 
repeats  itself  over  and  over  on  the  water's  surface.  Figure  1  shows 
overhead  photographs  of  two  such  hexagonal  patterns.  These  waves 
are  easy  to  generate  in  the  laboratory,  and  they  seem  to  be  stable  to 
perturbations. 

•A  nonlinear  partial  differential  equation  due  to  Kadomtsev  & 
Petviashvili  (1970), 

9x(9tU  +  6  u3xU  +  9x^u)  +  3  By^u  =  0,  (KP) 

is  known  to  describe  approximately  the  evolution  of  waves  in 
shallow  water.  The  equation  is  completely  integrable,  and  Krichever 
(1977)  showed  that  it  admits  huge  families  of  exact  solutions, 
including  an  eight-parameter  family  of  solutions  (called  "genus  2 
solutions")  that  are  spatially  periodic  and  have  peimanent  form. 
Segur  &  Finkel  (1985)  conjectured  that  this  family  of  KP  solutions 
ought  to  describe  all  spatially  periodic  waves  of  permanent  form  in 
shallow  water  of  uniform  depth.  Our  experimental  results  support 
that  conjecture.  We  have  developed  an  algorithm  to  find  the  eight 
parameters  of  the  KP  solution  that  best  fits  given  wave-gauge 
measurements  for  a  spatially  periodic  wave  of  permanent  form.  For 
every  experiment  we  ran,  the  KP  solution  so  obtained  described  the 
measured  waves  with  reasonable  accuracy,  even  when  the  waves 
violated  some  of  the  hypotheses  underlying  KP  theory.  As  an 
example.  Figure  2  shows  the  contour  lines  of  a  KP  solution  of  genus  2 
superposed  on  the  overhead  photgraphs  of  an  observed  periodic 
waves  in  shallow  water.  As  a  second  example.  Figure  3  shows 
twenty  seconds  of  data  from  an  array  of  nine  wave  gauges, 
measuring  the  surface  displacements  of  a  spatially  periodic  wave  of 
permanent  form  as  it  swept  by  the  gauges.  Also  shown  in  Figure  3 
are  the  corresponding  surface  displacements  predicted  by  the  KP 
solution  of  genus  2  that  best  fits  these  data.  Both  figures  show  good 
agreement  everywhere,  but  not  perfect  agreement.  The  reader  can 
find  more  information  about  this  study  in  Appendix  A. 
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2.  Work  in  progress 

(i)  Obliquely  reflected  waves 

The  waves  shown  in  Figures  1-3  were  created  with  a  specially 
designed,  complicated,  segmented  paddle  at  the  Waterways 
Experiment  Station  in  Vicksburg,  MI  (see  Scheffner,  1988,  for  a 
description  of  the  experimental  facility  there).  The  complexity  of 
this  facility  makes  one  wonder  if  the  waves  require  so  much 
complexity  for  their  generation,  or  if  simpler  mechanisms  might  also 
produce  spatially  periodic  waves  of  permanent  form. 

Figure  4  suggests  another  mechanism  to  produce  these  waves. 
The  figure  is  a  mosaic  of  two  simultaneous,  overhead  photographs  of 
the  water's  surface,  illuminated  by  a  light  shining  towards  the 
wavemaker.  In  the  figure,  a  train  of  nonlinear  plane  waves  (i.e.,  a 
cnoidal  wavetrain)  entering  the  picture  from  the  upper  left,  and 
reflects  obliquely  off  a  straight  vertical  wall  on  the  far  right.  The 
incident  and  reflected  waves  interact  to  produce  a  two-dimensional 
wave-pattern  that  is  spatially  periodic  (over  a  limited  spatial  region) 
and  that  propagates  along  the  wall  with  little  apparent  change  of 
form.  Because  oblique  reflections  of  plane  waves  are  common,  such  a 
mechanism  would  suggest  that  stationary,  hexagonal  wavepatterns 
like  those  shown  in  Figures  1-3  might  also  be  common. 

With  this  mechanism,  the  incoming  cnoidal  wave  should 
determine  the  outgoing  hexagonal  wave,  and  one  should  be  able  to 
deduce  of  the  parameters  of  the  outgoing  wave  from  those  of  the 
incoming  wave.  One  would  also  like  an  estimate  of  how  much 
distance  down  the  wall  is  required  for  this  transformation  from 
(incoming)  cnoidal  waves  to  (outgoing)  hexagonal  waves. 

Detailed  analysis  of  the  wave-records  for  the  waves  shown  in 
Figure  4  indicate  that  the  hexagonal  waves  observed  in  the  figure  are 
slowly  evolving  as  they  propagate  along  the  wall.  [This  evolution  is 
not  obvious  in  the  figure,  but  it  can  be  seen  in  the  corresponding 
wave-records.]  Apparently  the  test  section  in  the  experiment  shown 
in  the  figure  is  not  long  enough  for  the  waves  to  reach  a  steady  state. 
Plans  to  repeat  these  experiments  in  a  different,  simpler  but  longer 
facility  are  now  under  way.  (See  (iv)  below.) 

(ii)  Analysis  of  ocean  data 

The  KP  equation  predicts  the  existence  of  spatially  periodic 
waves  of  permanent  form  in  shallow  water  of  uniform  depth,  and 
laboratory  experiments  show  reasonably  good  agreement  with  these 
predictions.  However,  neither  of  these  statements  necessarily 
implies  that  waves  like  these  occur  in  real  oceans.  For  this,  one 
needs  actual  oceanic  measurements.  Fortunately,  oceanographers  at 
the  U.S.  Army  Engineer  Field  Research  Facility  (FRF)  in  Duck  N.C.  have 

3 


spent  several  years  measuring,  recording  and  cataloguing  the  waves 
incident  on  their  relatively  flat  beach.  Among  their  data  are 
measurements  from  an  array  of  eight  gauges,  located  about  915  m. 
offshore,  in  about  8  m.  of  water. 

In  July,  1994,  we  went  to  FRF  in  search  of  wave-records 
showing  spatially  periodic  waves  of  permanent  form.  Analysis  of  the 
data  is  still  in  progress,  but  the  first  comparison  is  shown  in  Figure  5. 
[We  are  grateful  to  Dr.  Chuck  Long  of  FRF  for  allowing  us  to  use  his 
data  for  these  comparisons.]  The  figure  shows  the  surface 
displacements  measured  at  eight  gauge  locations,  with  about  100 
seconds  of  data  at  each  gauge.  Also  shown  are  the  predictions  of  the 
best  KP  solution  for  these  data.  The  agreement  is  not  perfect,  and  it 
is  worse  than  that  for  the  laboratory  experiments.  In  particular,  the 
KP  solution  seems  to  get  the  phase  information  about  right,  but  it 
consistently  underestimates  the  peak  amplitudes. 

Even  so,  we  find  the  results  in  Figure  5  quite  encouraging. 

These  are  actual  ocean  waves,  measured  in  the  field,  with  no 
averaging  of  wave-records.  Figure  5  shows  that  a  single  KP  solution 
of  genus  2  simultaneously  approximates  the  data  from  all  eight 
gauges  with  reasonable  accuracy.  The  comparison  suggests  that  a 
comparatively  simple  model  might  explain  some  of  the  complex 
wave  patterns  observed  nearshore.  This  work  is  continuing. 

(iii)  Waves  of  genus  3 

All  of  our  work  in  shallow  water  to  date  has  involved  spatially 
periodic  waves  of  permanent  form,  approximated  by  KP  solutions  of 
genus  2.  However,  the  KP  equation  is  integrable,  so  it  also  admits 
more  complicated  solutions.  In  particular,  Krichever  (1977)  showed 
that  the  equation  admits  explicit,  quasiperiodic  solutions  of  any 
(integer)  genus.  A  solution  of  genus  N  represents  a  quasiperiodic 
function  with  N  independent  phases.  For  N  ^  3,  these  solutions  are 
typically  unsteady  in  every  coordinate  system,  so  they  ought  to 
provide  interesting  predictions  about  the  time-dependent  dynamics 
of  nonlinearly  interacting  waves  in  shallow  water. 

In  practice,  a  KP  solution  of  genus  N  is  defined  by  an  N-fold, 
nested  Fourier  series,  and  just  calculating  the  solutions  efficiently  is 
nontrivial.  In  fact,  no  general  algorithm  is  known  that  enables  one  to 
truncate  the  N-fold  series  to  within  a  predetermined  bound  on  the 
truncation  error.  A  more  complicated  problem  is  to  infer  the 
parameters  of  the  KP  solution  from  given  measurements  of  the  wave. 
Bobenko  &  Bordag  (1989)  have  developed  a  numerically  efficient 
procedure  to  calculate  these  series  for  certain  ranges  of  the 
parameters.  Whether  their  procedure  works  over  the  entire  range  of 
parameters  is  much  less  clear.  We  have  developed  an  alternative 


procedure  for  waves  of  genus  3,  which  is  based  on  earlier  work  of 
Dubrovin  (1981),  and  which  we  are  now  implementing  numerically. 
The  disadvantage  of  our  procedure  is  that  it  is  restricted  to  waves  of 
genus  1,  2  or  3,  while  that  of  Bobenko  &  Bordag  is  not  so  restricted. 

Its  advantages  are  that  it  provides  rigorous  bounds  on  truncation 
errors,  and  that  it  computes  the  KP  solution  with  as  few  terms  as 
possible.  A  more  comprehensive  comparison  of  the  two  procedures 
awaits  the  completion  of  the  numerical  code  for  our  procedure. 

(iv)  Experimental  work 

An  8 -axis  wavemaker  control  and  power  system  has  been 
configured,  purchased,  and  integrated.  The  system  comprises: 

(1)  8  brushless  DC  motor/drive  packages,  each  capable  of  a 
continuous  torque  of  14.2  oz.in.  at  speeds  up  to  175  rev/s;  and 

(2)  a  Programmable  Multi-Axis  Controller  (PMAC)  from  Delta  Tau 
Systems.  The  controller  is  a  state-of-the-art  system  that  exploits  DSP 
(digital  signal  processor)  technology  to  obtain  a  significant 
improvement  in  motion-control  over  traditional  servo-controllers. 

Commanded  motion  profiles  for  each  axis  can  be  calculated 
within  the  PMAC,  which  has  a  sophisticated  real-time  multi-tasking 
operating  system,  or  they  can  be  loaded  into  its  dual-ported  memory 
along  a  VME  bus.  The  controller  uses  encoder  feedback  from  each 
motor  to  obtain  the  instantaneous  position,  velocity  and  acceleration 
of  each  axis's  motion.  This  information  is  used  in  on-the-fly 
computations  to  command  the  motor  for  the  next  profile  move. 

The  8 -axis  motion  control  will  be  used  to  drive  eight  piston- 
type  wavemakers  in  an  existing  laboratory  basin.  (The  mechanical 
wavemaker  apparatus  is  still  under  development.)  The  need  for  this 
new  facility  is  twofold.  First,  the  large  facility  at  the  U.S.  Army 
Coastal  Engineering  Research  Center  (Vicksburg,  MS)  is  no  longer 
available.  Second,  the  high  degree  of  control  enables  us  to  scale  the 
experiments  down  in  size  so  that  more  flexibility  is  permitted  in 
planning  and  changing  experiments  (at  much  lower  costs).  Moreover, 
the  smaller  scale  will  enable  us  to  use  an  existing  high-speed  imaging 
system  that  measures  the  time  evolution  of  the  water  surface 
topography.  These  measurements  are  particularly  important  for 
investigating  2-D  water-wave  theories  such  as  the  KP  equation. 

The  first  set  of  experiments  planned  for  the  new  wavemaker 
system  involve  oblique  reflections  of  plane  waves,  as  discussed 
above.  The  second  set  of  experiments  planned  involve  investigations 
of  KP  solutions  of  genus  3,  also  discussed  above.  Because  KP 
solutions  of  genus  3  do  not  have  permanent  form,  the  surface 
imaging  measurements  will  be  particularly  useful  for  comparing  with 
KP  theory. 
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B,  Deep  water 
1 .  Pi^lished  work 

Numerical  simulations  provide  much  of  the  current  knowledge 
of  the  dynamics  of  waves  in  deep  water.  Consequently,  we  have 
focussed  our  attention  to  date  on  the  reliability  of  numerical 
simulations.  Our  work  has  appeared  in  five  published  papers: 

(i)  Ablowitz,  M.J.,  Schober,  C.M.  and  Herbst,  B.  Phys.  Rev.  Lett.,  71, 
pp.  2683-2686,  1993. 

(ii)  Ablowitz,  M.J.  and  Schober,  C.M.  "Effective  Chaos  in  the 
Nonlinear  Schrodinger  Equation"  in  Contemporary  Mathematics 
172,  Eds.  P.Kloeden  and  J.Palmer,  AMS  1994. 

(iii)  Ablowitz,  M.J.  and  Schober,  C.M.  "Hamiltonian  Integrators  for 
the  Nonlinear  Schroedinger  Equation"in  Computational  Physics. 
Ed.  J.  Potvin,  1994. 

(iv)  Ablowitz,  M.J.  and  Schober,  C.M.  "Numerical  Stochasticity, 

Hamiltonian  Integrators  and  the  Nonlinear  Schroedinger 
Equation",  in  the  Three  Dimensional  Dynamical  Systems 
Workshop.  Ed.  Dr.Kandrup,  N.  Y.  Academy  of  Sciences,  1993. 

(v)  Ablowitz,  M.J.  and  C.M.  Schober,  "Homoclinic  Manifolds  and 
Numerical  Chaos  in  the  Nonlinear  Schroedinger  Equation",  in 
Math,  and  Computers  in  Simulation  on  Nonlinear  Wavelike 
Eq*ns.  Eds.  Vichnevelsky,  Taha,  Newel.  Elsevier  Science,  1993. 

The  first  of  these  is  included  in  this  report  as  Appendix  B.  Our  main 
results  are  the  following. 

•  The  (focussing)  nonlinear  Schrodinger  equation  in  one  spatial 
dimension, 

iatA  +  ax2A  +  2  IA|2a  =  0,  (NLS) 

arises  in  several  physical  contexts,  including  as  an  approximate 
evolution  equation  for  the  complex  envelope  of  one-dimensional, 
nearly  monochromatic  waves  in  deep  water  (Zakharov,  1968).  The 
equation  is  known  to  be  integrable  (Zakharov  &  Shabat,  1972),  .so  its 
solutions  are  not  chaotic. 

•  Even  so,  for  certain  parameter  ranges  of  the  intial  data,  small 
errors  on  the  order  of  roundoff  (i.e.,  10*^^)  grow  rapidly  in  numerical 
simulations  and  saturate  at  values  comparable  to  the  main  wave. 
Because  NLS  is  integrable,  it  admits  several  constants  of  motion,  and 
these  constants  are  preserved  to  a  very  high  accuracy  during  the 
simulations.  Thus,  the  growth  of  these  numerical  errors  cannot  be 
attributed  to  standard  problems  like  inadequate  numerical  resolution 
or  poor  time-integration  strategies. 
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•  Instead,  the  numerical  errors  grow  because  of  a  phase  instability  in 
the  numerical  simulations.  The  instability  is  associated  with 
homoclinic  solutions  nearby  the  computed  solution.  Similar 
instabilities  occur  in  other  integrable  equations  that  admit  homoclinic 
solutions  (e.g.,  focussing  mKdV,  sine-Gordon),  but  not  in  integrable 
equations  that  admit  no  homoclinic  solutions  (e.g.,  defocussing  NLS, 
defocussing  mKdV,  KdV,  KP). 

•  The  instability  leads  to  chaotic  behaviour  in  the  numerical 
simulations,  so  one  observes  "numerically  induced  chaos".  In  other 
words,  even  if  NLS  itself  is  not  chaotic,  a  perturbed  NLS  equation  is 
chaotic,  and  tiny  numerical  errors  produce  enough  chaos  to 
completely  change  the  character  of  the  solution.  This  example 
demonstrates  the  danger  of  inferring  that  a  Hamiltonian  dynamical 
system  is  chaotic  on  the  basis  of  numerical  simulations  of  the  system. 

•  To  relate  these  results  to  problems  in  water  waves,  we  note  that 
the  sign  difference  between  focussing  and  defocussing  NLS  is  the  sign 
difference  that  determines  whether  or  not  the  Benjamin-Feir 
instability  is  operative.  Thus,  we  would  expect  numerically-induced 
chaos  in  numerical  simulations  of  waves  in  deep  water,  but  not  for 
waves  in  shallow  water. 

2.  Work  in  progress 

For  focussing  NLS,  by  changing  a  parameter  (M)  in  the  initial 
data  we  can  excite  an  arbitrary  number  of  linearly  unstable  modes. 
When  M=2-3  we  have  observed  that  at  intermediate  values  of  the 
mesh  size,  the  standard  numerical  schemes  break  down  and  develop 
a  serious  spatial  and  temporal  chaotic  response  (Herbst  &  Ablowitz, 
1989,  Ablowitz  &  Herbst,  1993).  However,  we  have  developed  an 
algorithm  that  is  based  upon  the  Inverse  Scattering  Transform  and 
that  we  call  the  integrable  discrete  NLS  equation  or  IDNLS,  which  is 
fast  and  effective.  In  this  parameter  regime  the  IDNLS  scheme  does 
not  show  any  signs  of  the  joint  spatial-temporal  chaos  that  exists  in 
the  standard  schemes.  When  the  mesh  is  sufficiently  refined,  the 
standard  schemes  eventually  converge  to  the  same  values  as  those 
obtained  via  the  IDNLS  scheme,  although  they  require  considerably 
finer  mesh  size. 

When  the  number  of  unstable  modes  is  larger,  e.g.  M=6-8,  we 
have  found  that  the  numerical  simulations  over  moderate  time 
intervals  are  extremely  sensitive.  It  is  difficult  to  compute  accurate 
solutions  even  though  we  use  an  extremely  refined  mesh  and  our 
numerical  algorithm  accurately  preserves  the  underlying  constants 
of  motion  (Ablowitz,  Herbst  &  Schober,  1993).  We  find  that 
numerical  instabilities  and  chaos  are  excited  even  by  tiny  errors  on 
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the  order  of  roundoff.  The  reason  is  that  the  initial  data  in  this 
parameter  regime  is  extremely  close  to  the  homoclinic  manifolds  of 
the  NLS  equation.  In  order  to  show  this,  one  needs  to  analyze 
carefully  the  inverse  spectral  theory  associated  with  the  NLS 
equation  (Ablowitz  &  Schober,  1993).  Perturbations  in  the 
neighborhood  of  homoclinic  manifolds  allows  frequent  homoclinic 
crossings  and  the  ensuing  chaotic  dynamics.  On  the  other  hand, 
when  we  consider  the  defocusing  NLS  equation,  theory  tells  us  that 
there  is  no  linear  instability  of  the  basic  wavetrain  and  there  are  no 
associated  homoclinic  manifolds.  In  this  case  our  numerical 
algorithms  work  well  and  there  is  no  indication  of  a  chaotic  response. 
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Figure  Captions 

Figure  1.  Mosaic  of  two  overhead  photographs,  showing  surface 
patterns  of  waves  in  shallow  water.  Experiment:  (a)  KP400203- 
300303;  (b)  KP220203- 110403.  Each  of  these  wave  patterns  has  a 
basic  hexagonal  template;  one  such  template  is  shown  in  Figure  la. 
Figure  2.  For  experiment  KP220203- 110403,  contour  lines  of  the 
best  KP  solution  are  superimposed  on  the  photograph  shown  in 
Figure  lb.  The  direction  of  propagation  of  the  entire  pattern  is  also 
shown. 

Figure  3.  Detailed  comparison  of  measured  wave  records  with  the 
best  KP  solution,  for  20  seconds  of  data  from  experiment  KP220203- 
110403. 

Figure  4.  Mosaic  of  two  overhead  photographs,  showing  a  cnoidal 
wavetrain  from  the  upper  left  impinging  obliquely  on  a  vertical  wall 
at  the  far  right.  The  incident  and  reflected  waves  interact  to  produce 
a  hexagonal  wave  pattern  that  has  nearly  permanent  form  as  it 
propagates  along  the  wall  towards  the  bottom  of  the  photo. 

Figure  5.  This  figure  shows  105  seconds  of  data,  measured  on  23 
September,  1990  at  8  wave  gauges  situated  in  8  m.  of  water  offshore 
of  FRF  in  Duck,  N.C.  Also  shown  superposed  on  these  data  are  the 
predictions  of  a  KP  solution  of  genus  2.  The  KP  solution  was  chosen 
to  fit  the  data,  and  this  solution  is  compared  to  the  data  from  all  8 
gauges  simulatneously. 
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Abstract: 

We  demonstrate  experimentally  the  existence  of  a  family  of 
gravity-induced,  finite-amplitude  water  waves  that  propagate 
practically  without  change  of  form  in  shallow  water  of  uniform 
depth.  The  surface  patterns  of  these  waves  are  genuinely  two- 
dimensional,  and  periodic.  The  basic  template  of  a  wave  is 
hexagonal,  but  it  need  not  be  symmetric  about  the  direction  of 
propagation,  as  required  in  our  previous  studies  (c.g.,  Hammack  et 
aU  1989).  Like  the  symmetric  waves  in  earlier  studies,  the 
asymmetric  waves  studied  here  are  easy  to  generate,  they  seem  to 
be  stable  to  perturbations,  and  their  amplitudes  need  not  be  small. 
The  Kadomtsev-Petviashvili  (KP)  equation  is  known  to  describe 
approximately  the  evolution  of  waves  in  shallow  water,  and  an  eight- 
parameter  family  of  exact  solutions  of  this  equation  ought  to  describe 
almost  all  spatially  periodic  waves  of  permanent  form.  We  present 
an  algorithm  to  obtain  the  eight  parameters  from  wave-gauge 
measurements.  The  resulting  KP  solutions  are  observed  to  describe 
the  measured  waves  with  reasonable  accuracy,  even  outside  the 
putative  range  of  validity  of  the  KP  model. 
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1.  Introduction: 

An  earlier  paper  by  Hammack,  Scheffner  &  Segur  (1989, 
hereinafter  called  Part  1)  reported  an  experimental  study  of  gravity- 
induced  waves  that  propagate  practically  without  change  of  form  in 
shallow  water  of  uniform  depth.  Those  waves  had  finite  amplitudes, 
and  their  surface  patterns  were  periodic  in  two  spatial  dimensions 
and  in  time.  The  spatial  pattern  of  the  waves  was  a  hexagon  that 
was  symmetric  about  the  direction  of  wave  propagation.  The  waves 
were  predicted  with  reasonable  accuracy  by  a  family  of  exact 
solutions  of  an  equation  due  to  Kadomtsev  and  Petviashvili  (1970): 

(ft  +  6  ffx  +  fxxx)x  +  3fyy  =  0,  (KP) 

where  subscripts  denote  partial  derivatives.  In  a  second  set  of 
experiments  (Hammack,  Scheffner  &  Segur,  1991),  more  effort  was 
made  to  reduce  depth-variations  in  the  laboratory  basin,  and  the 
theoretical-experimental  agreement  improved  further.  Both  of  these 
studies  were  preceded  by  Peregrine  (1985),  who  also  observed 
hexagon-like  waves  of  permanent  form  but  who  argued  against  that 
interpretation. 

Within  the  KP  equation,  almost  all  real-valued,  spatially 
periodic  waves  of  permanent  form  are  characterized  by  a  family  of 
exact  solutions  with  eight  free  parameters:  six  dynamic  parameters 
associated  with  wavelengths  and  wave  amplitudes,  plus  two 

arbitrary  phases  (c/.  Segur  &  Finkel,  1985).  Typically,  these  waves 
have  asymmetric  surface  patterns.  The  two  experimental  papers 
cited  above  tested  only  a  subset  of  these  waves,  called  symmetric 
waves,  which  have  only  (3  +  2)  parameters.  Figure  1  shows 
photographs  of  both  symmetric  and  asymmetric  water  waves. 

In  this  paper,  we  describe  an  experimental  study  of  two- 
dimensional,  spatially  periodic,  asymmetric  waves  that  propagate 
practically  without  change  of  form  in  shallow  water  of  uniform 

depth.  [In  the  terminology  used  in  this  paper,  surface  waves  are 
either  "one-dimensional"  or  "two-dimensional",  depending  on  the 
nature  of  their  surface  patterns.  The  velocity  fields  associated  with 
the  waves  are  one  dimension  higher.]  Our  three  main  results  for 
asymmetric  waves  are  consistent  with  our  earlier  results  for 

symmetric  waves: 

1)  In  shallow  water  of  uniform  depth,  there  exist  spatially  periodic 
waves  with  finite  amplitudes  and  fully  two-dimensional  surface 

patterns.  These  waves  are  easy  to  generate,  they  propagate  with 
nearly  permanent  form,  and  they  seem  to  be  stable  with  respect  to 
perturbations  in  initial  conditions  and  to  weak  external  forcing  (e.g.. 
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due  to  viscosity  or  to  variations  in  bottom  depth).  The  amplitudes  of 
the  waves  need  not  be  small. 

2)  The  basic  template  of  the  two-dimensional  wave  pattern  is  a  six- 
sided  figure  (i.e.,  a  hexagon),  with  a  broad,  flat  trough  surrounded  by 
six  narrow  wave  crests  (e.g.,  see  Figure  1).  Wave  crests  opposite 
each  other  in  this  hexagon  are  parallel,  they  are  equal  in  length  and 
in  wave  amplitude,  and  they  can  be  identified  with  each  other  in 
order  to  reproduce  the  periodic  pattern. 

3)  Over  the  range  of  parameter  values  that  we  tested,  these  water 
waves  are  predicted  with  reasonable  accuracy  by  the  full,  eight- 
parameter  family  of  KP  solutions  of  genus  2.  This  range  extends  well 
beyond  the  putative  range  of  validity  of  the  KP  equation  as  an 
approximate  model  of  waves  in  shallow  water  (i.e.,  beyond  the  range 
of  nearly  one-dimensional  waves  of  small  amplitude  in  shallow 
water).  In  this  sense,  the  KP  model  provides  a  convenient  means  to 
learn  about  hexagonal  waves  of  permanent  form,  but  their 
mathematical  existence  and  stability  seem  to  be  independent  of  the 
KP  model. 

Our  primary  objective  in  this  study  is  to  establish  experimentally 
the  existence  and  properties  of  asymmetric  hexagonal  waves  of 
permanent  form  in  shallow  water  of  uniform  depth.  To  this  end,  we 
present  the  results  of  fifteen  experiments,  which  exhibit  some  of  the 
variety  of  waves  of  this  form.  Our  secondary  objective  is  to 
determine  how  accurately  the  full,  eight-parameter  family  of  KP 
solutions  of  genus  2  describe  the  measured  waves.  To  this  end,  we 
present  an  algorithm  that  uses  measurements  from  an  array  of 
wave-gauges  to  obtain  the  eight  parameters  of  the  KP  solution  that 
"best"  fit  the  Pleasured  waves.  By  design,  the  algorithm  uses  only 
wave-gauge  measurements  (ignoring  the  known  data  used  to 
program  the  wavemaker)  so  that  it  can  be  applied  to  gauge-array 
data  in  more  complex  wave  environments  in  the  laboratory  and  in 
the  field.  Fifteen  experiments  are  far  too  few  to  determine  the 
parameter  range  for  KP  theory  to  be  applicable,  and  our  experiments 
were  not  chosen  with  that  objective  in  mind.  In  fact,  most  of  our 
experiments  have  parameters  that  lie  outside  the  putative  range  of 
validity  of  KP  theory.  (Some  of  our  waves  have  near-beaking  heights 
and  strongly  two-dimensional  surface  patterns.)  We  reiterate  our 
belief  that  the  KP  model  is  convenient  but  not  essential  to  describe 
periodic  waves  of  permanent  form  in  shallow  water. 

A  summary  of  the  contents  of  this  paper  is  as  follows.  In  §2,  we 
review  the  main  ideas  needed  to  compare  KP  solutions  of  genus  2 
with  experimental  data  on  waves  in  shallow  water.  We  describe 
experimental  procedures  and  results  in  §3.  The  procedures  are 
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similar  to  those  used  in  Part  1,  but  the  experimental  data  here  are 
more  complicated.  The  main  theoretical  problem  addressed  in  this 
paper  is  the  following.  Given  wave-gauge  data  from  a  particular 
experiment,  find  the  eight  parameters  of  the  KP  solution  of  genus  2 
that  best  fits  these  data.  We  present  an  algorithm  to  solve  this 
problem  in  an  Appendix  to  the  paper,  and  we  use  the  algorithm  in  §4 
to  compare  KP  theory  with  the  results  of  fifteen  different 
experiments.  In  a  sense,  our  procedure  can  be  viewed  as  a 
complicated  form  of  surface-fitting  (i.e.,  two-dimensional  curve¬ 
fitting),  with  no  demonstrated  predictive  power.  However,  we 
recorded  more  data  than  we  used  to  find  the  KP  solution,  and  in  §4 
we  also  test  the  accuracy  of  some  concrete  predictions  of  the  theory. 
The  hexagonal  waves  studied  in  this  paper  and  in  part  1  seem  to  be 
the  simplest  nontrivial  waves  of  permanent  form  that  are  periodic 
and  two-dimensional  in  shallow  water.  This  identification  suggests 
that  they  ought  to  occur  frequently  in  natural  oceanic  settings.  In  §5, 
we  present  some  observations  of  oceanic  waves  similar  to  those 
discussed  in  this  paper. 


2.  Review  of  KP  theory 
The  KP  equation. 


(ft  +  6  ffx  +  fxxx)x  +  3fyy  =  0,  (KP) 

describes  approximately  the  slow  evolution  of  gravity-induced 
waves  on  water  of  uniform  depth  when  the  waves  are  assumed; 
(i)  to  have  horizontal  lengthscales  much  longer  than  the  fluid  depth 
(this  assumption  is  also  called  "shallow  water");  (ii)  to  propagate 
primarily  in  one  direction  (the  x-direction),  v.ith  only  weak 
variations  in  the  transverse  y-direction;  and  (iii)  to  have  small-to- 
moderate  amplitudes.  The  KP  equation  is  not  well-posed  as  it  stands 
(f  H  0  and  f  =  t  both  solve  KP,  and  they  coincide  at  t  =  0),  and  one 
must  interpret  {3x'M  in  order  to  integrate  (KP)  in  time.  All  of  the  KP 
solutions  discussed  in  this  paper  satisfy  a  constraint: 

|f(x,y,t)dx  =  0, 

and  this  constraint  eliminates  the  ambiguity  of  interpreting  • 

Ablowitz  and  Villarroel  (1991)  studied  this  and  related  constraints. 
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A  detailed  derivation  of  the  KP  equation  can  be  found  in  Segur 
&  Finkel  (1985),  or  elsewhere.  Briefly,  one  begins  with  Stokes' 
(1847)  equations  for  surface  waves  on  an  inviscid,  incompressible 
fluid  resting  on  a  horizontal  bed  under  a  constant  gravitational  force, 
g.  Imposing  the  three  assumptions  above  yields  at  leading  order  the 
linear,  one-dimensional  wave  equation,  whose  solution  consists  of 
left-going  and  right-going  waves,  each  travelling  with  speed  igh, 
where  h  is  the  mean  depth  of  the  fluid.  Each  set  of  waves,  of  course, 
has  permanent  form  at  this  order  (i.e.,  for  short  times).  Over  longer 
times,  in  a  coordinate  system  moving  with  (say)  the  right-going 
waves,  one  observes  that  the  right-going  waves  evolve  slowly,  due  to 
the  three  small  effects  that  were  neglected  at  leading  order.  The  last 
three  terms  in  the  KP  equation  represent  the  effect  on  this  slow 
evolution  of  (i)  weak  nonlinearity,  (ii)  weak  dispersion,  and  (iii)  weak 
two-dimensi  onality . 

Let  {X,  Y,  Z)  represent  spatial  coordinates  in  a  fixed  laboratory 
frame  with  Z  vertical,  let  T  be  time,  let  'n(X,  Y,  T)  measure  the 
elevation  of  the  fluid  surface  above  its  mean  level,  and  let  e  «  1  be  a 
formal  small  parameter  related  to  the  assumptions  above.  The 
scaled  variables  in  (KP)  can  be  related  to  the  laboratory  coordinates 
as  follows: 


x  =  ii(K-ighT)/h, 
t  =  ^f^X/6h, 


y  =  eY/h, 

f  =  3€Ti/2h  +  0(e2). 


(2) 


This  change  of  coordinates  can  be  performed  in  more  than  one  way; 
the  choice  used  here  is  slightly  unconventional,  and  it  differs  from 
that  in  Part  1:  ordinarily  one  sets  t  =  e^^V^T/6h.  Either  choice  leads 
to  a  KP  equation,  but  the  solution  must  be  interpreted  slightly 
differently  in  the  two  cases.  The  current  choice  is  more  natural  for 
our  experimental  configuration,  but  the  main  motivation  for  this 
choice  is  that  it  simplifies  the  surface-fitting  algorithm  presented  in 
the  Appendix. 

As  noted  above,  the  KP  equation  is  an  approximate  model  of 
water  waves  for  e  «  1.  In  order  to  compare  with  definite 
experiments,  we  set  e  =  1  in  this  paper.  To  justify  this  prescription 
(e=l),  we  note  that  the  KP  equation  is  invariant  under  the  scaling: 

X  ->  ax,  y  a^y,  t  ->  a^t,  f  a'^f.  (3) 
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Applying  this  scaling  to  (2)  without  the  O(e^)  terms,  with  =  e, 

amounts  to  setting  e  =  1  there.  As  a  consequence,  some  free 
parameters  in  the  solutions  described  below  must  be  small  for  the 
validity  of  the  approximation. 

Krichever  (1977)  showed  that  the  KP  equation  admits  a  huge 
family  of  exact  solutions  in  the  form 

f(x,  y.  t)  =  2  ax^ln  e,  (4) 

where  0  is  a  Riemann  theta  function  of  genus  N.  Here  N  identifies 
the  number  of  independent  phases  in  the  solution,  and  f  is  a  quasi- 
periodic  function  of  these  N  phases.  (Recall  that  "quasi-periodic" 
means  that  f  is  periodic  in  each  of  the  N  phases,  if  the  other  N-1 
phases  are  held  fixed.)  In  the  simplest  case,  N  =  1  and  (4) 
reproduces  the  familiar  cnoidal  wave  (e.^.,  Wiegel,  1960).  The 
solutions  of  interest  in  this  paper  have  N  =  2,  and  6  is  given  by  a 

double  Fourier  series: 

••  •• 

6(^i.<I>2;  B)  =  X  X  +  2mnbX  +  n2(bxV  d))  +  i(nM|>i  +  n^l,  (5) 

11= 

where  the  phases  and  <>2  are  given  by 

<|»j  =  PjX  +  Vjy +  o)jt  +  <Dj,  (j  =  l,  2)  (6) 

( 1 ,  02}  are  arbitrary  phase  constants,  and  the  real-valued 
parameters  (b,  X,  d),  which  define  the  elements  of  the  Riemann 
matrix  (see  Segur  and  Finkel,  1985),  satisfy  the  following  constraints: 

-oo<b<0,  0<X<1,  -00  < d  ^  Kl-X^).  (7) 

Actually,  the  constraint  on  X  can  be  tightened  to  [0  <  X  ^  1/2].  Even 
so,  we  use  (7)  because  it  allows  us  to  assume  that  pi/p2^  0,  which 
simplifies  the  algorithm  in  the  Appendix.  For  solutions  of  genus  2, 
formal  asymptotic  validity  of  the  KP  model  of  water  waves  requires 
the  limit  {pj  0,  with  (vj/pj^)  and  (oj/pj^)  finite).  However,  we  will 
demonstrate  reasonable  agreement  between  theory  and  experiment 
even  for  values  of  {pi,  P2}  that  are  not  very  small. 

The  eight  free  parameters  of  a  real-valued  KP  solution  of  genus 
2  are  {b,  X,  d,  pi,  p2,  V2;  d>i,  d»2}.  As  noted  by  Dubrovin  (1981)  and  by 
Segur  &  Finkel  (1985),  the  parameters  can  be  selected  as  follows: 
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(i)  choose  {b,  X,  d}  to  satisfy  (7); 

(ii)  choose  1x2  as  an  arbitrary  positive  number,  and  v  2  as  an  arbitrary 
real  number; 

(iii)  choose  (pi/p2)  so  that  a  sixth-order  polynomial  in  ()Xi/|i2)  is 
positive  (the  polynomial  can  be  found  in  Dubrovin  (1981)  or  in  Segur 
&  Finkel  (1985));  then  the  square  root  of  the  polynomial  yields  vi; 

(iv)  once  vi  is  known,  then  (Oi  and  <02  can  be  found  explicitly  in  terms 
of  {Pi,li2}; 

(v)  choose  {<]>i,d>2}  as  arbitrary  real  numbers,  but  since  they  are 
phase  constants,  it  is  sufficient  to  permit  only 

0  ^  Oj  <  2jc.  0  =  1.  2) 

Every  such  choice  of  parameters  yields  a  KP  solution  of  genus  2 
that  is  real-valued,  that  is  periodic  in  two  spatial  directions  and  in 
time,  and  that  propagates  without  change  in  form.  The  velocity  of 
propagation  (in  KP  coordinates)  is  given  by 

dx  _  Vl0)2  -  V2<0l  dy  _  Pi0)2  -  M-20>1 

dt  ■  P1V2  -  P2V1  ’  dt  "  p2Vi  •  fiiV2 '  (8 ) 

We  obtain  KP  solutions  using  a  FORTRAN  code  named  DELTA.  The 
code  is  available  to  interested  readers  through  an  anonymous  ftp 
site;  see  Ap{)endix  A  for  details. 

Symmetric  solutions  of  genus  2  are  obtained  by  imposing  three 
additional  restrictions  on  these  parameters: 

d  =  b(l-X2),  (Pi/H2)  =  1.  V2  =  -V,.  (9) 

It  follows  from  these  that  0D2  =  ti>i>  and  that  the  KP  solution  is 
invariant  under  the  transformation  {y  -»  -y,  x  x,  t  -»  t).  Among  the 
KP  solutions  of  genus  2,  symmetric  waves  are  special  in  at  least  three 
ways: 

(i)  a  symmetric  wave  propagates  purely  in  the  x-direction,  because 
dy/dt  =  0  in  (8); 

(ii)  a  symmetric  wave  is  specified  by  only  (3+2)  free  parameters, 
instead  of  by  (6+2),  because  of  (9);  and 

(iii)  a  symmetric  wave  is  periodic  in  the  x-  and  y-directions. 
Symmetric  solutions  were  the  main  focus  of  the  work  in  Part  1,  but 
not  here. 
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3.  Experimental  program  <t  typical  results 

The  experimental  program  was  similar  to  that  described  in 
Part  1  and  identical  to  that  described  briefly  in  Hammack,  Scheffner 
&  Segur  (1991).  Experiments  were  performed  at  the  Coastal 
Engineering  Research  Center  (CERC),  U.S.  Army  Engineer  Waterways 
Experiment  Station,  Vicksburg  MS.  We  used  a  basin,  whose  planform 
is  shown  schem  cally  in  Figure  2,  that  comprised  a  uniform-depth 
section  12.5Sm  long  and  26.52m  wide,  and  a  section  with  a  gently 
sloping  (1:30)  beach  that  absorbed  most  of  the  incident  wave  energy. 
In  the  uniform-depth  section,  the  still  water  depth  was  20cm  +3 mm. 

A  segmented  wavemaker  with  58  piston-type  paddles,  each 
45.7  cm  wide,  spanned  the  wide  wall  opposite  the  beach.  The 

wavemaker  was  described  in  detail  by  Scheffner  (1988),  and  the 

wavemaker  motion  was  described  in  Part  1.  Briefly,  the  wavemaker 
was  programmed  to  generate  a  superposition  of  two  cnoidal 
wavetrains  with  different  directions  of  propagation.  In  the  current 

experiments,  we  allowed  the  wavelengths  and/or  amplitudes  of  the 

two  cnoidal  waves  to  differ,  so  that  the  resulting  two-dimensional 
wave  patterns  would  be  asymmetric.  The  nomenclature  of  our 
experiments,  e.g.  KPaabbcc-AABBCC,  is  based  on  programming 
parameters  for  the  wavemaker.  The  aa  and  AA  represent  the  phase 
lags  in  degrees  between  adjacent  paddles  of  the  segmented 
wavemaker  for  each  of  the  cnoidal  wavetrains;  the  bb  and  BB  are  the 
wavelengths  in  meters  for  each  of  the  cnoidal  wavetrains;  and  the  cc 
and  CC  are  crest  amplitudes  (i.e.,  maximum  elevations  above  the 
mean  water  level)  in  ems  for  each  of  the  cnoidal  wavetrains. 
Ordinarily,  the  wavelengths  and  amplitudes  of  the  two-dimensional 
wave  patterns  that  resulted  from  these  inputs  differed  slightly  from 
these  values.  A  1  Os-ramp  was  applied  to  the  start  of  each  command 
signal  driving  the  wavemaker  to  protect  the  mechanical  system  from 
sudden  starts.  Hence,  there  is  a  10s  transient  period  at  each  site  in 
the  basin  before  the  programmed  waves  arrived. 

Our  primary  means  of  wave  measurement  was  a  linear 
array  of  nine  wave  gauges  that  was  oarallel  to,  and  8m  from,  the 
wavemaker  astride  the  basin  centerline.  (See  Figure  2.)  The  nine 
capacitance-type  gauges,  spaced  at  2m  intervals,  measured  the 
elevation  of  the  water  surface.  A  second  array  of  similar  gauges  was 
placed  normal  to  the  wavemaker  along  the  basin  centerline.  It 
contained  eight  gauges  (one  shared  with  the  first  array)  in  the 
uniform-depth  section;  seven  spaced  at  0.5m  intervals  beginning 
6.5m  from  the  wavemaker,  and  one  2m  from  the  wavemaker.  (This 
array  also  contained  two  gauges  in  the  sloping-beach  section;  data 
from  these  gauges  are  not  reported  here.)  The  wave  gauges  were 
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calibrated  and  recording  began  in  a  quiescent  basin  immediately 
prior  to  starting  the  wavemaker.  Continuous-time  signals  for  lOOs 
were  obtained  from  ail  18  gauges;  these  signals  were  low-pass 
filtered  using  an  8-pole  butterworth  filter  with  a  lOHz  cutoff,  and 
then  digitized  to  produce  25Hz  discrete-time  data.  The  algorithm  to 
fit  KP  solutions  to  wave  data,  described  in  detail  in  Appendix  B,  is 
based  solely  on  data  from  the  nine-gauge  array;  i.e.,  we  make  no  use 
of  the  wavemaker  parameters,  which  are  given  in  the  experiment 
label.  As  stated  earlier,  the  algorithm  was  designed  in  this  way  so 
that  it  could  also  be  used  for  waves  of  unknown  origin.  Data  from 
the  second  array  are  used  in  §4  to  test  predictions  of  the  fitted 
solutions. 

Wavefields  were  also  measured  qualitatively  using  two 
Hasselblad  cameras,  placed  6m  apart  astride  the  basin  centerline  and 
7m  above  the  nine-gauge  array.  The  cameras  took  simultaneous 
overlapping  pictures  of  the  water  surface  in  a  darkened  room  using 
four  880  W-s  strobe  lights  located  on  the  beach  and  pointing  toward 
the  wavemaker.  (The  wave  gauges  were  removed  from  the  basin 
during  photography.)  The  two  photographs  were  then  combined  to 
form  a  mosaic,  like  those  in  Figure  1.  All  photographs  were  printed 
so  as  to  maintain  the  same  horizontal  scale,  which  is  indicated  by  the 
(+)  marks  painted  on  the  basin  floor  at  Im  intervals.  The  small,  but 
measurable,  differences  in  spacing  between  (+)  marks  apparent  in 
Figure  1  resulted  from  optical  distortions  by  *he  waves,  which  acted 
like  moving  lenses.  This  distortion  prevented  us  from  using  the 
overhead  photographs  to  measure  spatial  wavelengths 
quantitatively. 

Fifteen  experiments  on  asymmetric  waves  were  conducted  in 
the  set  under  discussion.  In  each  experiment,  a  two-dimensional, 
spatially  periodic  wave  pattern  was  created  at  one  end  of  the  tank. 
The  entire  wave  pattern  propagated  across  the  tank  with  nearly 
permanent  form,  it  was  recorded  as  it  swept  by  the  array  of  gauges, 
and  then  it  was  mostly  absorbed  on  the  sloping  beach. 

Figure  1  shows  clear  photographic  evidence  that  two- 
dimensional,  approximately  periodic  waves  exist  in  shallow  water. 
Photographs  of  three  different  wave  patterns,  corresponding  to  three 
experiments,  are  shown.  Each  wave  pattern  propagates  from  top  to 
bottom  in  the  photos;  the  front  of  each  wave  crest  appears  bright, 
while  the  back  appears  dark.  Some  of  the  bright-dark  boundaries  in 
Figure  la  are  particularly  sharp,  indicating  that  these  wave  crests 
are  particularly  steep,  and  close  to  breaking.  Close  inspection  of 
Figure  la  also  shows  that  on  the  front  faces  of  these  steep  waves  are 
capillary  waves  that  appear  to  be  radiating  from  the  crests,  as 
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expected  for  near-breaking  waves.  Thus,  even  though  the  KP 
equation  is  asymptotically  valid  only  for  waves  of  small  amplitude. 
Figure  la  clearly  demonstrates  that  two-dimensional,  approximately 
periodic  waves  exist  even  at  large  amplitudes.  The  existence  and 
stability  of  these  waves  seem  to  persist  well  beyond  the  putative 
range  of  validity  of  the  KP  model. 

The  wave  pattern  in  Figure  la  is  symmetric,  like  those  discussed 
in  Part  1,  and  the  pattern  propagates  directly  downward  in  the 
photograph,  i.e.  in  the  x-direction.  The  wave  pattern  in  Figure  lb  is 
slightly  asymmetric,  and  that  in  Figure  Ic  is  strongly  asymmetric. 
Asymmetric  waves  propagate  in  directions  oblique  to  the  x-direction. 
A  feature  common  to  all  of  the  waves  in  all  of  our  experiments, 
including  those  in  Figure  1,  is  the  basic  hexagonal  shape  of  the  wave 
pattern.  The  basic  template  has  a  broad,  flat  trough  surrounded  by 
six  relatively  narrow  wave  crests.  One  such  hexagon  is  drawn  in 
Figure  lb.  The  wave  crests  opposite  each  other  in  this  hexagon  are 
parallel,  they  are  equal  in  length  and  in  wave  amplitude,  and  they 
can  be  identified  with  each  other  in  order  to  reproduce  the  periodic 
pattern,  i.e.,  one  can  tile  the  plane  with  copies  of  this  basic  hexagon. 

Figure  3  shows  typical  data  obtained  from  the  nine-gauge 

array  for  three  experiments  on  asymmetric  waves,  including  two 
that  were  photographed  in  Figures  Ib,c.  (Note  that  the  waves  shown 
in  Figures  lb,c  represent  about  3  seconds  of  data  in  Figures  3a,b.) 
The  100s  records  show  an  initially  quiescent  water  level  and  the 
subsequent  arrival  of  the  generated  wave  pattern,  whose  first  10s 
exhibit  a  transient  caused  by  the  ramped  start  of  the  wavemaker. 

Data  taken  during  these  intervals  are  excluded  in  subsequent 

analysis.  (We  note  also  that  some  wave  reflection  undoubtedly 
occurs  on  the  1:30  sloping  beach,  leading  to  small  standing-wave 
components  in  the  wave  records.  No  attempt  has  been  made  to 
remove  reflected  waves  from  the  data.)  At  first  glance,  data  at  each 
gauge  site  may  appear  to  be  periodic,  especially  in  Figure  3b; 

however,  careful  examination  shows  that  the  periodicity  is  not  exact, 
which  is  characteristic  of  quasi-periodic  waves. 

Figure  4  shows  the  amplitude-frequency  spectrum  of  each 
of  the  wave  records  shown  in  Figure  3b.  (The  spectral  amplitude  is 
the  modulus  of  the  fast  Fourier  transform,  or  FFT,  of  the  data.  For 
convenience,  we  will  refer  to  these  amplitudes  as  FFTs.]  To  save 
space,  we  omit  the  spectra  for  the  other  experiments  shown  in  Figure 
3,  but  these  spectra  are  qualitatively  similar  to  those  shown  here. 
Figure  4  shows  clearly  the  characteristic  signatures  of  two- 
dimensional,  nonlinear,  two-phase  waves:  each  of  the  three  FFTs 
shows  two  dominant  frequencies,  plus  their  harmonics  and  some 
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sum-  and  difference-frequencies.  As  one  moves  from  the  FFT  of  one 
wave  record  to  another  within  each  experiment,  the  amount  of 
energy  (i.e.,  square  of  the  spectral  amplitude)  at  each  frequency 
changes,  but  the  total  energy  is  always  distributed  among  the  same 
two  frequencies,  plus  their  harmonics.  In  this  way,  the  FFTs  for  each 
experiment  show  that  all  of  the  gauges  measured  one-dimensional 
slices  through  the  same  two-dimensional  wave  pattern.  The 

photographs  in  Figures  lb,c  show  clearly  the  spatial  pattern  of  two- 
dimensional,  two-phase  waves,  so  one  expects  temporal  FFTs  like 
those  in  Figure  4.  However,  based  solely  on  the  FFTs  of  wave-gauge 
data  we  could  infer  the  following  information  directly: 

(i)  the  wave  pattern  has  two  phases  (or  possibly  more  than 
two  in  degenerate  cases,  but  with  only  two  independent 
frequencies  among  the  phases); 

(ii)  the  wave  pattern  is  nonlinear,  because  energies  at 

harmonics  as  well  as  at  the  sum-  and  difference- 

frequencies  are  significant. 

The  fact  that  the  same  two  frequencies  are  measured  at  every  gauge, 
with  different  energy  levels,  suggests  that  the  wave  pattern  might  be 
two-dimensional,  but  it  does  not  establish  the  two-dimensionality 
conclusively.  If  the  wave  pattern  is  genuinely  two-dimensional  with 
exactly  two  phases,  then  necessarily  it  is  stationary  in  some 

uniformly  translating  coordinate  system;  i.e.,  it  must  be  a  wave  of 
permanent  form. 

Figure  4  also  lists  the  numerical  values  of  the  two  fundamental 
frequencies  measured  in  the  experiment.  The  two-dimensional  wave 
patterns  shown  in  Figure  1  are  periodic,  but  the  one-dimensional 

wave  records  in  Figure  3  are  not  periodic  unless  these  two 

frequencies  are  rationally  related.  Since  the  frequencies  listed  in 

Figure  4  were  obtained  numerically,  they  are  necessarily  rationally 

related.  However,  for  experiment  KP220203-1 10403,  the  period 
implied  by  the  two  frequencies  exceeds  600  seconds,  so  this  exact 
periodicity  is  irrelevant  for  the  100  seconds  of  recorded  data,  which 
are  effectively  quasi-periodic. 

We  note  finally  that  quasi-periodic  functions  also  arise  in  the 
theory  of  dynamical  systems  (c.g.,  Guckenheimer  &  Holmes,  1983), 
where  quasi-periodic  behaviour  often  suggests  incipient  instability 
and  chaos  (Piexoto,  1962).  No  such  instability  is  suggested  from  the 
quasi-periodicity  seen  in  Figure  3.  These  data  are  quasi-periodic  only 
because  the  data  are  taken  from  a  doubly  periodic  function,  along  a 
line  at  an  irrational  angle  from  a  direction  of  periodicity.  The  water 
waves  being  measured  seem  to  be  remarkably  stable. 


4.  Comparison  of  theory  and  experiment 

In  Appendix  B,  we  present  an  algorithm  to  find  the  KP  solution  of 
genus  2  that  best  Bts  experimental  data  like  those  shown  in  Figure  3. 
Here  is  a  brief  summary  of  the  algorithm.  Recall  from  §2  that  a  KP 
solution  of  genus  2  is  specified  by  eight  real-valued  parameters. 
Four  of  these  (pi,  P2.  vi,  vi)  are  wavenumbers  that  determine  the 
overall  spatial  structure  of  the  wave  pattern.  Two  of  them  are  easy 
to  find:  (pi,  P2)  are  proportional  to  the  two  frequencies  identified  in 
Figure  4.  Once  the  spatial  structure  is  known,  then  one  varies  (b,  d, 
X.)  in  order  to  match  (as  closely  as  possible)  the  observed  maximum 
and  minimum  amplitudes,  and  the  known  spatial  structure.  Then 
one  determines  the  two  phase  constants  (Oi,<I>2)  by  minimizing  the 
rms  (or  L2)  error  for  the  entire  wave  pattern.  Finally,  with  good 
approximations  to  the  eight  KP  parameters  determined  in  this  way, 
one  selects  final  values  for  all  eight  parameters  by  minimizing  the 
rms-error  for  the  wave  pattern  again.  (Thus,  the  solution  found 
in  this  way  is  "best"  in  the  sense  of  least  rms-error.) 

This  algorithm  was  used  to  find  the  best  KP  solution  for  each  of 
fifteen  experiments  with  asymmetric  wave  patterns.  We  begin  our 
discussion  of  results  with  the  question:  How  accurately  does  the  best 
KP  solution  represent  the  data  used  to  select  that  solution?  For 
conciseness,  we  concentrate  on  the  three  experiments  shown  in 
Figure  3. 

a)  Surface-fitting 

For  the  experiment  shown  in  Figures  lb  and  3a,  Figure  5  shows  a 
comparison  between  the  measured  data  with  the  best  KP  solution  at 
the  same  locations,  for  20  seconds  of  data,  and  for  all  of  the  gauges  in 
the  nine-gauge  array.  (Most  of  our  comparisons  begin  30  seconds 
into  the  recorded  data  in  order  to  exclude  the  quiescent  and 
transient  intervals,  and  last  for  20  seconds.  The  agreement  between 
theory  and  experiment  is  largely  independent  of  when  we  start  the 
comparison  and  how  long  we  compare.)  The  agreement  between  KP 
theory  and  experimental  data  in  Figure  5  is  perfect  nowhere,  but  it  is 
fairly  good  everywhere.  The  KP  solution  captures  the  overall 
structure  of  the  wave  pattern;  it  describes  the  phase  information, 
such  as  zero-crossings,  particularly  well.  Note  especially  the 
amplitudes  of  the  waves  in  Figure  5.  From  Figure  3,  the  maximum 
crest-to-trough  height  (H)  for  this  wave  pattern  is  11.23  cm.,  which 
corresponds  to  a  wave  height-to-depth  ratio  (H/h)  of  0.56.  Peregrine 
(1983)  discusses  criteria  for  wave-breaking  in  shallow  water  (c.g., 
H/h  =  0.77  at  breaking).  The  wave  pattern  shown  here  is  not 
breaking,  but  it  is  not  far  from  the  range  of  breaking  waves; 
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certainly  the  wave  amplitudes  are  not  especially  small.  This  result  is 
similar  to  that  found  in  Part  1  for  symmetric  waves:  KP  theory 
agrees  well  with  experimental  measurements  well  outside  its 
putative  range  of  validity. 

Figure  6  shows  another  comparison  of  theory  and  experiment, 
for  the  wave  pattern  shown  in  Figures  Ic  and  3b.  This  pattern  is 
strongly  asymmetric,  but  again  the  best  KP  solution  gives  good 
agreement  everywhere.  Figure  3b  shows  that  these  waves  are 
almost  periodic,  with  an  almost-period  of  about  20  seconds.  Figure  6 
shows  that  the  KP  solution  has  essentially  the  same  behaviour.  [In 
fact,  the  comparison-time  of  20  seconds  was  chosen  to  contain  the 
almost-period  of  these  waves.]  Because  this  experiment  was  also 
photographed  (in  Figure  Ic),  we  can  compare  the  KP  solution  with 
the  photograph  as  well.  In  Figure  7,  contour  lines  of  the  best  KP 
solution  are  superimposed  on  the  photograph.  This  comparison 
contains  different  information  than  that  in  Figure  6:  X  is  held  fixed 
while  T  and  Y  vary  in  Figure  6,  whereas  T  is  held  fixed  while  X  and  Y 
vary  in  Figure  7.  Moreover,  the  photograph  was  not  used  to  find  the 
KP  solution,  so  Figure  7  is  not  surface-Htting.  Even  so.  Figure  7  leads 
to  a  conclusion  similar  to  that  obtained  from  Figures  S  and  6:  the  best 
KP  solution  describes  the  observed  wave  pattern  with  reasonable 
accuracy  everywhere. 

The  last  detailed  comparison,  for  KP3(X)1 02-30*803,  has  been 
included  because  the  wave  records  from  gauges  3  and  7  showed  a 
slow  increase  in  amplitude  over  the  course  of  this  experiment,  as 
seen  in  Figure  3c.  While  it  is  clear  that  the  data  in  Figure  3b  are 
almost  periodic,  this  is  much  less  clear  in  Figure  3c,  so  this 
experiment  provides  a  strong  test  of  the  theory.  Figure  8a  shows  a 
comparison  of  the  measured  wave  records  with  the  best  KP  solution, 
for  20  seconds  of  data  starting  at  30  seconds.  These  waves  have 
smaller  amplitudes  than  those  in  other  experiments  in  this  series,  so 
the  waves  are  more  nearly  sinusoidal,  and  the  agreement  between 
theory  and  experiment  is  slightly  better  than  average.  Note  that  the 
wave  amplitudes  for  gauges  3  and  7  increase,  even  in  this  20-second 
interval.  Note  also  that  for  gauge  3,  the  KP  solution  captures  this 
behaviour. 

Figure  8b  shows  a  second  comparison  for  the  same  experiment, 
KP300 102-30*803,  starting  at  70  seconds.  It  is  important  to  note 
that  the  KP  solution  used  here  is  the  same  KP  solution  as  that  in 
Figure  8a,  except  that  the  two  phase  variables,  Oi  and  d>2,  were  re¬ 
optimized  for  the  new  data,  which  were  measured  during  a  different 
time  interval.  In  other  words.  Figure  8b  is  not  a  surface-fit;  it  is  a 
prediction  of  the  theory,  based  on  the  surface-fit  shown  in  Figure  8a. 
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Figure  3c  shows  that  the  data  measured  at  20  seconds  differ 
considerably  from  those  measured  at  70  seconds,  but  the  same  KP 
solution  describes  both  data  sets  with  reasonable  accuracy.  This 
good  agreement  implies  that  the  data  in  Figure  3c  are  almost 
periodic,  but  the  almost-period  is  so  much  longer  than  the  80 
seconds  of  data  that  this  feature  is  not  apparent  in  the  measured 
data.  We  mention  that  we  also  applied  the  optimization  algorithm  to 
the  data  in  Figure  8b,  to  find  the  l^st  KP  solution  for  those  data.  The 
KP  solution  found  was  nearly  identical  with  that  obtained  from  the 
data  in  Figure  8a,  except  for  changes  in  <I>i  and  (I>2- 
b)  Analysis  of  errors 

Table  1  summarizes  the  comparisons  between  KP  theory  and 
experiments  for  fifteen  experiments  involving  asymmetric  waves. 
The  parameters  of  the  best  KP  solutions  are  listed,  as  well  as  the 
maximum  and  minimum  measured  wave  amplitudes,  and  two 
measures  of  error,  which  we  now  define. 

For  a  given  experiment,  let  u(x,  y,  t)  denote  the  measured  wave 
amplitude  at  location  (x,  y,  t),  normalized  as  in  (2),  and  let  f(x,  y,  t) 
denote  the  value  of  the  best  KP  solution  at  the  same  location.  One 
measure  of  their  discrepancy  is  o,  where 

9 

X  I  -  u(x,y,t)l2dx 

q2  -  I  Jo _ 

9 

Xj  [u(x,y,t)l2dx 

the  sum  is  taken  over  the  wave  records  from  the  nine-gauge  array, 
and  L  is  the  length  of  the  measured  wave  records.  (We  generally 
used  20  seconds  of  data,  but  varying  the  length  of  the  data  string  did 
not  significantly  affect  o.)  A  second  measure  of  error  is  fmaxl^max^ 
the  ratio  of  maximum  wave  amplitudes  according  to  theory  and 
experiment.  (Step  4  of  the  optimization  algorithm  effectively 
minimizes  (I-/max/«max)^»  but  then  the  last  step  of  the  algorithm 
usually  forces  it  away  from  its  minimal  value,  leaving  f  maxi  Umax  as  a 
nontrivial  measure  of  error.) 

According  to  Table  1,  the  KP  model  along  with  the  optimization 
algorithm  consistently  underpredicted  the  maximum  wave 
amplitude.  This  error  was  more  serious  in  some  experiments  than  in 
others,  but  the  model  predicted  a  maximum  amplitude  that  was  too 
small  in  every  experiment.  This  consistent  underprediction  was 
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apparently  caused  by  using  the  norm  in  (10).  In  Part  1  we  used  a 
different  norm,  and  matched  the  maximum  amplitudes  well.  Beyond 
this  simple  observation,  we  found  no  striking  trends  in  the  data.  For 
all  fifteen  experiments,  o  stayed  within  a  fairly  narrow  range: 

0.235  ^  a  ^  0.402,  with  a  mean  value  of  E(a)  =  0.302.  (For 
comparison,  note  that  the  trivial  theory,  in  which  f(x,  y,  t)  s  0,  gives 
0=1  in  (10).)  We  do  not  suggest  from  this  relatively  narrow  range  of 
values  that  the  theory  never  breaks  down,  but  only  that  our 
experiments  did  not  demonstrate  the  breakdown. 

As  mentioned  in  §2,  the  KP  model  of  water  waves  is  based 
on  three  assumptions:  (i)  moderate  wave  amplitudes;  (ii)  shallow 
water;  and  (iii)  nearly  one-dimensional  wave  propagation.  An 
important  question  is  how  a  depends  on  each  of  these  three  effects. 
Consider  first  the  wave  amplitudes,  which  we  characterize  by  the 
maximum  crest-to-trough  height,  divided  by  the  undisturbed  water 
depth,  which  was  20  cm  in  these  experiments.  (This  choice  is  based 
on  the  columns  entitled  Vmax  Vmin  Table  1,  and  is  somewhat 
arbitrary.  We  also  tried  using  Tlmax*  similar  results.  Using  a, 

instead  of  or  some  other  measure  of  error,  is  also  arbitrary.)  For 
the  experiment,  let 

Aj  :=  (llmax,  j  ■  'Hmin,  jl/l*.  (11) 

From  Table  1,  0.341  ^  Aj  ^  0.561,  with  a  mean  value,  E(A)  =  0.492. 
The  correlation  coefficient  between  A  and  a  in  our  fifteen 
experiments  is 


p(A,  a)  := 


X  [Aj-E(A)]*[aj-E(a)] 

15  15 

(X  [Aj-E(A)]2*X  [Oj-E(a)]2}i/2 


0.39. 


(12) 


This  positive  correlation  shows  that  KP  theory  becomes  less  accurate 
as  wave  amplitudes  increase,  as  expected. 

It  remains  to  test  the  assumptions  of  shallow  water,  and  of 
nearly  one-dimensional  waves.  Unfortunately,  our  data  set  (based 
on  fifteen  experiments)  was  too  small  to  carry  out  these  tests 
conclusively. 
c)  Predictions 

Knowing  that  the  KP  model  provides  reasonably  accurate 
descriptions  of  the  measured  wavefields,  we  now  turn  to  other 
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comparisons  between  theory  and  experiments.  For  each  experiment, 
the  best  KP  solution  is  now  fixed,  with  no  remaining  free  parameters; 
hence,  the  comparisons  that  follow  involve  predictions  of  the  theory, 
rather  than  just  surface-fitting. 

A  fundamental  assertion  in  this  paper  is  that  the  waves  under 
study  propagate  with  nearly  permanent  form.  Yet  data  from  the 
nine-gauge  array  cannot  test  this  assertion,  because  all  nine  gauges 
are  the  same  distance  from  the  wavemaker.  Limited  insight  into 
possible  wave  evolution  is  provided  by  the  second  gauge  array, 
which  was  placed  perpendicular  to  the  first  array  (see  Figure  2).  We 
show  next  that  to  within  our  experimental  error,  the  waves  did  not 
evolve  significantly  over  the  7.5  m  span  of  the  second  array.  Again, 
for  simplicity,  we  concentrate  on  the  three  experiments  shown  in 
Figure  3. 

For  each  experiment,  the  second  gauge-array  provides  eight 
wave  records  like  those  in  Figure  3.  Figure  9  shows  the  temporal  FFT 
of  each  of  the  seven  new  wave  records  for  experiment  KP220203- 
110403,  featured  in  Figures  Ic,  3b,  and  4.  (Note  from  Figure  2  that 
gauge  S  was  included  in  both  arrays.  In  presenting  the  data  from 
the  second  array,  we  omit  gauge  S  because  its  data  were  already 
presented.  We  also  omit  the  data  from  the  other  experiments,  which 
data  were  qualitatively  similar  to  those  shown  in  Figure  9.)  Note 
that  every  gauge  in  the  second  array  recorded  the  same  two 
dominant  frequencies,  and  therefore  the  same  values  for  pi  and  \i2’, 
in  other  words,  pi  and  p2  showed  no  measurable  evolution  as  the 
waves  propagated  away  from  the  wavemaker. 

Figure  10  shows  the  actual  wave  records  from  the  seven  new 
gauges,  for  experiment  KP220203- 110403.  The  KP  solution  that  was 
already  selected,  based  on  data  from  the  nine-gauge  array,  predicts 
the  data  that  each  of  these  seven  gauges  should  measure;  these 
predictions  are  also  shown  in  Figure  10.  We  emphasize  that  no  free 
parameters  were  available  in  the  comparisons  in  Figure  10;  these  are 
predictions  rather  than  surface-fits.  The  overall  error  (ox)  for  the 
seven-gauge  array  is  comparable  to  that  for  the  nine-gauge  array  (o), 
even  though  data  from  the  latter  array  were  used  to  select  the  KP 
solution.  For  all  fifteen  experiments  on  asymmetric  waves.  Table  2 
lists  an  overall  error  (a)  for  the  nine-gauge  array  used  to  select  the 
KP  solution,  and  another  (ox)  for  the  seven-gauge  array,  whose 
measurements  were  predicted  by  that  KP  solution.  Within  each 
experiment,  these  errors  are  comparable;  usually  ax  is  slightly  larger 
than  o,  but  for  experiment  KP400203_300303,  the  error  for  the 
seven-gauge  array  (ax)  is  actually  smaller  than  that  for  the  nine- 


gauge  array  (o).  This  comparison  supports  the  claim  that  the 
measured  waves  did  not  evolve  significantly  over  the  length  of  our 
test  section. 

Approximate  wavelengths  in  the  X-direction  for  a  KP  solution 
can  be  obtained  from  {pi,p2)>  using  Li ->  2nh/pi.  For  experiment 
KP220203_ 110403,  this  gives  X-wavelengths  of  about  2  m  and  4  m, 
consistent  with  the  cnoidal  wavelengths  input  to  the  wavemaker. 
Note  that  gauge  10  is  only  2  m.  from  the  wavemaker,  so  Figure  10 
shows  that  the  wave  being  measured  achieved  its  apparently 
permanent  form  within  one  wavelength  of  the  wavemaker,  and  then 
retained  that  form  at  least  for  the  length  of  the  test  section  of  the 
tank  (7.5  m).  In  all  fifteen  experiments,  we  always  found  that  the 
wave  patterns  achieved  their  apparently  permanent  form  within  2  m 
of  the  wavemaker,  and  then  retained  that  form  over  the  length  of 
our  test  section.  We  found  no  clear  evidence  of  wave  evolution. 

Next,  we  consider  the  (constant)  velocity  of  propagation  of  the 
waves.  This  velocity  cannot  be  inferred  from  the  photographs  in 
Figure  1.  The  KP  model  predicts  it,  in  (8),  and  the  velocity  vector 
shown  in  Figure  7  was  found  in  this  way.  However,  the  velocity  also 
can  be  obtained  directly  from  the  wave  records,  without  using  (8),  as 
follows.  For  a  particular  experiment,  Hnd  the  KP  solution  that  best 
fits  the  data  measured  during  a  20  second  interval  starting  at  time 
(T).  If  the  wave  has  permanent  form,  then  at  a  slightly  later  time 
(T  +  AT),  this  wave  pattern  will  have  translated  in  space,  and  what 

had  been  observed  at  location  (X  -  AX,  Y  -  AY)  in  the  fixed 

(laboratory)  frame  now  will  be  observed  at  (X,  Y),  for  some  {AX,  AY). 
In  order  to  measure  (AX,  AY),  use  (2)  to  replace  (6)  by 

<|>j  =  pj(x  -  AX/h)  +  Vj(y  -  AY/h)  +  a)j(t  -  AX/6h)  +  d>j,  (j  =  1,2)  (13) 

and  then  minimize  at  the  new  time  (i.e.,  over  20  seconds  starting 

at  T  +  AT)  by  varying  {AX,  AY),  holding  fixed  all  the  other  parameters 
of  the  solution.  In  this  way,  one  obtains  empirically  {AX,  AY)  for  a 
given  AT.  For  a  wave  of  permanent  form,  (AX/AT)  and  (AY/AT)  do 
not  depend  on  AT,  and  it  follows  from  (2)  that 

AX  AX  ^  ^ 

Because  the  KP  solution  is  spatially  periodic,  {AX,  AY)  are  ambiguous 
in  the  sense  that  one  can  add  to  them  an  integer  number  of  spatial 
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periods  of  the  wave  pattern.  One  eliminates  this  discrete  ambiguity 
by  requiring  that  AX  ->  0,  AY  -»  0  as  AT  0;  then  {AX,  AY)  are 
unambiguous  for  AT  small  enough.  On  the  other  hand,  small 
measurement  errors  become  relatively  important  if  AT  is  too  small. 
Balancing  these  competing  effects,  we  chose  AT  =  1.2  seconds  for 
these  experiments. 

Table  2  lists  the  two  components  of  the  velocity  of 
propagation  of  the  wave  pattern,  according  to  KP  theory  (8)  and 
according  to  direct  measurement  (14),  for  all  fifteen  experiments. 
The  agreement  is  remarkably  good;  in  many  cases,  the  two  velocities 
agree  to  within  a  few  percent.  The  agreement  becomes  even  more 
impressive  when  one  recalls  that  these  velocities  are  themselves  the 
small  deviations  from  the  basic  wave  speed,  Vgh,  which  is  inherent  in 
the  (moving)  KP  coordinate  system. 
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5.  Oceanic  observations  of  hexagonal  waves  in  shallow  water 

The  central  point  of  this  paper  is  the  following.  In  shallow 
water  of  uniform  depth,  there  exist  two-dimensional  periodic  waves 
with  finite  amplitudes  and  nearly  permanent  form.  The  spatial 
pattern  of  these  waves  is  hexagonal,  as  shown  in  Figure  1.  The 
waves  are  easy  to  generate  in  the  laboratory,  and  they  are 
apparently  stable. 

KP  solutions  of  genus  2  describe  these  waves  with  reasonable 
accuracy  over  a  range  of  parameters  that  is  larger  than  the  putative 
range  of  validity  of  the  KP  model.  The  robustness  of  the  KP  model  is 
a  pleasant  surprise,  but  we  suspect  that  the  existence  of  these  waves 
is  independent  of  the  validity  of  the  KP  model.  Rather,  we  conjecture 
that  the  hexagonal  wave  pattern  observed  in  Figure  1  is  the  basic 
qualitative  pattern  for  periodic  waves  of  permanent  form  in  shallow 
water.  More  precisely,  we  conjecture  that  if  one  seeks  genuinely 
two-dimensional  wave  patterns  of  finite  amplitude  that  are  periodic 
and  that  propagate  with  permanent  form  in  shallow  water  of 
uniform  depth,  then  one  will  almost  certainly  find  hexagonal  wave 
patterns  with  flat  troughs  and  narrow  crests,  like  those  in  Figure  1. 
If  this  conjecture  is  correct,  then  one  would  expect  to  find  waves  like 
these  in  oceanic  settings.  In  this  section,  we  present  some  oceanic 
observations  suggesting  that  wave  patterns  like  these  actually  occur. 
Remark:  Some  apparent  counter-examples  to  this  conjecture  are 
actually  limiting  cases  of  it.  We  mention  specifically:  (i)  cnoidal 
waves,  in  which  the  wave  period  in  one  direction  becomes  infinitely 
long;  (ii)  waves  of  infinitesimal  amplitude,  for  which  the  hexagon 
degenerates  to  a  parallelogram,  and  the  waves  become  sinusoidal; 
and  (iii)  the  oblique  interaction  of  two  solitons,  which  is  not  periodic, 
but  is  another  limiting  case  of  hexagonal  waves. 

The  first  example,  shown  in  Figure  11,  is  an  aerial  photograph 
of  the  coastal  zone  taken  near  Jones  Inlet  on  Long  Island  near  New 
York  City.  At  the  top  of  the  photograph  one  can  identify  a  road,  some 
buildings,  the  beach,  and  the  surf  zone  in  which  waves  are  breaking. 
Beyond  the  surf  zone,  in  the  bottom  half  of  the  photograph,  one  sees 
clearly  that  the  surface  wave  patterns  are  two-dimensional  and 
approximately  periodic,  like  the  waves  in  Figure  1,  and  that  many  of 
the  wave  patterns  are  hexagonal. 

The  second  example,  shown  in  Figure  12,  is  a  now  well-known 
photograph  taken  by  Terry  Toedtcmcicr  of  two  waves  interacting  in 
shallow  water  off  the  Oregon  coast.  This  photograph  is  often  shown  to 
illustrate  the  interaction  of  two  solitons,  with  the  required  phase 
shift  (i.e.  a  spatial  shift  of  each  wave  crest  as  a  result  of  the  nonlinear 
interaction  of  the  two  waves).  However,  each  of  the  most  prominant 
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wave  crests  in  Figure  12  is  actually  part  of  a  train  of  periodic  waves. 
The  next  crest  in  each  train  is  discernible  in  the  figure,  but  the 
strongest  evidence  of  the  periodicity  of  the  wavefield  is  that  Mr. 
Toedtemeier  told  us  that  the  waves  were  periodic.  He  also  noted  that 
none  of  the  other  wave  interactions,  before  or  after  the  one  shown  in 
Figure  12,  was  nearly  as  dramatic  as  the  one  shown.  Thus,  these 
waves  were  only  approximately  periodic,  but  they  were  certainly 
two-dimensional,  and  they  exhibited  flat  troughs  and  sharp  crests.  In 
our  model,  the  two  long  crests  and  the  one  short  crest  shown  in  the 
figure  should  be  viewed  as  the  sides  of  two  adjacent  hexagons. 

The  last  example,  shown  in  Figure  13,  is  an  aerial  photograph  of 
waves  off  the  outer  banks  of  North  Carolina  taken  by  Carl  Miller 
during  a  major  storm  (the  "Hallowe’en  storm  of  1991").  The  waves  in 
this  photograph  were  enormous — wave  gauges  located  in  10m  deep 
water  measured  waves  with  heights  exceeding  Sm  and  periods  of 
20s.  If  we  take  the  wave  speed  to  be  Vih  as  a  rough  approximation, 
then  a  wave  period  of  20s  corresponds  to  a  wavelength  of  about 
200m.  (Assigning  a  width  of  5  m.  to  the  road  visible  along  the  coast 
and  then  scaling  also  suggests  wavelengths  on  the  order  of  100-200 
m.)  The  photograph  itself  suggests  that  virtually  every  wave  crest 
was  breaking,  and  thai  the  beach  on  the  right  of  the  photo  was 
completely  flooded.  These  waves  were  apparently  far  outside  the 
range  of  validity  of  the  KP  model.  Moreover,  these  waves 
undoubtedly  did  not  have  permanent  form,  so  that  KP  solutions  of 
genus  2  would  be  too  simple  to  describe  them.  Nevertheless,  the 
waves  in  Figure  13  exhibit  some  of  the  same  features  that  we  have 
now  seen  repeatedly:  two-dimensional  wave  patterns  with  flat 
troughs  and  sharp  crests,  and  with  many  of  the  crests  forming 
hexagonal  wave  patterns. 

The  photographs  in  Figures  11,  12  and  13  suggest  that 

approximately  periodic,  hexagonal  wave  patterns  occur  in  the  ocean. 
Without  a  systematic  study  of  oceanic  observations  of  waves  in 
shallow  water,  we  do  not  know  whether  these  waves  are  common  or 
uncommon.  We  suspect  that  they  are  more  common  than  is  usually 
believed. 
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Appendix  A:  To  obtain  the  program  DELTA 

The  KP  solutions  used  in  this  paper  were  computed  using  a 
FORTRAN  program  named  DELTA.  The  following  instructions  show 
how  to  use  the  UNIX  ftp  utility  to  transfer  the  file  "delta.f"  from  the 
anonymous  site  "ftp.colorado.edu:/cuboulder/appm".  The  indented 
lines  are  keyboard  commands,  with  <cr>  =  carriage  return. 

To  log  onto  the  ftp  site,  type 

ftp  ftp.colorado.edu  <cr> 

You  should  get  a  prompt  for  a  username.  Type 
anonymous  <cr> 

You  will  then  be  prompted  for  a  password.  Enter  your  complete 
email  address.  At  this  point,  you  should  be  logged  onto  the  ftp  site. 
You  need  to  go  to  the  subdirectory  /cuboulder/appm.  Type 
cd  cuboulder/appm  <cr> 

To  transfer  the  file,  type 
get  delta.f  <cr> 

To  quit  the  ftp  utility,  type 
quit  <cr> 
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Appendix  B:  An  algorithm  to  find  the  best  KP  solution  for  a  set  of 
wave-gauge  data 

A  KP  solution  of  genus  2  is  defined  by  (6+2)  parameters 
{b,  X,  d,  m.  112.  V2;  <bi,  O2),  from  which  the  other  parameters  of  the 
solution  (V],  0)1,  0)2)  can  be  deduced.  The  objective  of  this  Appendix 
is  to  provide  an  effective  algorithm  to  identify  a  "best"  set  of  KP 
parameters  for  a  given  set  of  wave  records,  like  those  in  Figure  3. 
The  algorithm  has  six  steps,  and  it  makes  explicit  use  of  the  fact  that 
the  time-variation  of  a  quasi-periodic  function,  like  that  of  (5)  and 
(6),  corresponds  to  straight-line  motion  on  a  two-dimensional  torus. 
Step  1:  Find  pj  and  P2- 

If  the  measured  wave  amplitudes  can  be  represented  by  a  KP 
solution  of  genus  2,  then  it  follows  from  §2  that  we  seek  a  function  of 
the  form 

f(x,  y,  t)  =  F(pix  +  viy  +  (oit,  it2x  +  V2y  +  t02t), 

where  F  is  periodic  in  each  argument.  Without  loss  of  generality,  we 
can  require  that  each  period  be  27c;  otherwise,  the  precise  form  of  F 
does  not  affect  step  1.  In  terms  of  laboratory  coordinates,  after  using 
(1)  with  e  =  1, 

f  =  F(pi(X  -  irihT)/h  +  ViY/h  +  0)iX/6h,  P2(X  -  VghT)/h  +  ViY/h  +  a>zXJ6h).  (15) 

Each  wave  gauge  records  data  as  T  varies,  with  X  and  Y  Hxed.  If  it 
records  a  function  in  the  form  (15),  then  the  FFT  of  the  wave  record 
shows  energy  at  two  frequencies,  Vi^Piand  V^P2,  plus  harmonics 
and  sum-  and  difference-frequencies.  In  other  words,  one  can  read 
Ipil  and  Ip2l  directly  from  FFTs  of  the  wave  records,  like  those  in 
Figure  4.  Because  of  (7),  one  can  always  require  p  1/112  ^  0-  Because 
0  is  even  [0(4>i,  <|>2)  =  -<t>2)].  one  can  always  take  pi  ^  0,  P2  ^  0- 

These  two  symmetries  allow  us  to  determine  pi  and  P2  from  FFTs. 
Remark:  The  advantage  of  the  choice  of  variables  in  (2)  is  that  it 
permits  direct  measurement  of  pi  and  P2-  With  the  usual  choice, 
t  =  e^^y^T/6h,  the  FFTs  yield  (pi  -©1/6)  and  (p2  *©2/6)  and  one  must 

solve  to  find  pi  and  P2-  We  have  tried  both  approaches;  the  final 
results  differ  only  slightly,  and  the  logic  of  the  current  approach  is 
simpler. 

Step  2:  Find  /v//  and  /V2/. 
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In  these  experiments  we  used  only  nine  gauges  at  different 
locations  in  Y,  so  a  Fourier  transform  in  Y  of  the  data  at  fixed  (X,  T) 
would  yield  only  crude  estimates  of  Ivil  and  Iv2l-  We  obtain  more 
precise  estimates  of  Ivil  and  Ivil  when  the  wave  height  is  a  quasi- 
periodic  function  of  two  variables,  as  follows. 

Let  F(<t>i,  <t>2)  be  continuous  and  quasi -periodic.  Along  the  line 
[<|)i  =  a<|>2  +  PI  in  the  ((^i,  <|)2)-plane,  F(4)i,  <j)2)  is  periodic  if  and  only  if  a 
is  rational.  If  a  is  irrational,  then  along  this  line  F(<|>i,  (|>2)  never 
repeats  exactly,  but  it  comes  arbitrarily  close  to  every  value  taken 
by  F(<|>i,  <|>2)  in  the  periodic  square.  Now  compare  this  situation  to 
that  in  Figures  1  and  3.  Figure  1  shows  a  periodic  function  of  two 
variables,  and  each  wave  record  shown  in  Figure  3  is  taken  along  a 
straight  line  defined  by  the  direction  of  travel  of  the  wave  pattern  as 
it  sweeps  past  the  gauges.  It  follows  that  if  the  data  in  each  record 
are  not  periodic,  then  a  sufficiently  long  string  of  data  from  a  single 
gauge  would  eventually  record  the  entire  two-dimensional  wave 
pattern,  to  any  desired  accuracy.  Moreover,  the  nine  gauges  sweep 
out  nine  parallel  line-segments  in  the  ((>1,  (|>2)-plane.  If  the  data  are 
not  periodic,  then  to  any  desired  accuracy,  one  can  think  of  these 
nine  line-segments  as  segments  on  the  same  (infinite)  line,  with 
different  starting  points. 

For  example,  in  Figure  3  b,  the  data  recorded  on  gauge  9  nearly 
coincides  with  that  recorded  on  gauge  3,  but  with  a  time-shift  of  36.8 
seconds.  This  can  be  seen  directly,  by  laying  a  record  from  gauge  9 
over  that  from  gauge  3  (but  shifted  by  36.8  seconds)  and  observing 
that  the  data  nearly  coincide.  More  quantitatively,  the  two  data  sets 
(with  the  first  20  seconds  of  each  data-set  deleted,  to  remove 
quiescent  and  transient  intervals)  have  a  correlation  coefficient  of 
p  =  0.966,  where 


P  = 


ni(t)  nj(t)  dt 


1/2 


(16) 


and  where  ni(t)  and  113(0  each  have  zero  mean.  Moreover,  using  the 
same  time-shift  (36.8  seconds),  the  data  from  gauge  8  correlate  with 
those  from  gauge  2  (with  p  =  0.966),  and  gauge  7  reproduces  the  data 
from  gauge  1  (with  p  =  0.975).  In  other  words,  with  an  error 
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corresponding  to  a  correlation  coefficient  of  about  0.97,  a  shift  in  Y  of 
12  m  (the  distance  between  correlated  gauges)  is  equivalent  to  a 
shift  in  T  of  36.8  seconds.  We  have  just  demonstrated  this  with 
three  pairs  of  gauges:  9->3,  8->2,  7-4  1,  By  this  reasoning,  had  there 
been  a  "virtual"  gauge  10  in  the  shore-parallel  array,  it  would  have 
measured  a  signal  (nearly)  equal  to  that  measured  by  gauge  4,  but 
shifted  by  36.8  seconds.  Moreover,  virtual  gauge  11  would  have 
measured  a  signal  corresponding  to  that  of  gauge  S,  and  so  on. 
Eventually,  when  we  return  to  gauge  9  (at  virtual  gauge  15),  we 
could  start  over  by  comparing  its  signal  to  that  of  gauge  3  again. 
This  entire  process  can  be  repeated  until  we  have  shifted  (in  time)  to 
the  end  of  the  data-records.  In  addition,  the  data  from  a  single 
gauge  can  be  matched  to  itself  in  a  second,  independent  way,  and 
between  these  two  shifts  we  can  extend  the  data  indefinitely. 
Referring  again  to  Figure  3  b,  note  that  the  data  from  any  one  gauge 
are  almost  periodic,  with  an  almost-period  of  41.84  seconds.  In 
particular,  in  Figure  3  b  the  correlation  coefficient  between  the  data, 
and  the  data  shifted  by  41.84  seconds,  averages  about  0.99  for 
gauges  4-9.  By  making  use  of  these  two  kinds  of  shifts,  one  can 
create  as  many  virtual  gauges  as  desired.  Each  shift  introduces  a 
small  error,  whose  size  is  related  to  (1  -  p).  In  principle,  one  should 
stop  adding  virtual  gauges  when  the  accumulation  of  these  errors 
exceeds  the  increased  precision  obtained  from  adding  more  gauges. 
Preliminary  tests  suggested  that  the  break-even  point  for  these  data 
occurs  at  about  200  gauges;  in  practice,  we  always  stopped  at  200 
gauges.  Including  the  9  original  gauges. 

Now  take  a  FFT  in  Y  at  fixed  (X,  T),  based  on  200  data;  the 
results  are  shown  in  Figure  14,  for  experiment  KP220203-1 10403. 
As  in  Figure  4,  the  energy  is  concentrated  in  two  dominant 
"frequencies",  plus  harmonics  and  sum-  and  difference-frequencies. 
Assuming  that  the  data  represent  a  function  of  the  form  (10),  then 
these  dominant  frequencies  are  located  at  Ivi/hl  and  Iv2/hl.  In  this 
way,  the  procedure  provides  direct  measurement  of  Ivil  and  IV2I.  This 
procedure  does  not  specify  which  of  the  two  values  should  be 
associated  with  pi,  and  which  with  P2-  Until  this  has  been  resolved 
(in  step  3),  let  us  call  the  two  values  Ival  and  Ivbl,  with  Ival  ^  Ivbl. 

Step  2a:  Symmetric  waves 

The  essential  requirements  for  step  2  are  that  the  data  be 
quasi-periodic,  but  not  strictly  periodic.  If  the  data  are  periodic  with 
a  sufficiently  short  period,  then  the  wave  records  from  two  different 
gauges  might  never  look  alike,  for  any  time-shift.  In  particular, 
symmetric  waves  produce  periodic  wave  records,  and  the  records 
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from  different  gauges  look  alike  only  if  two  of  the  gauges  happen  to 
be  an  integer  number  of  Y-wavelengths  apart.  For  symmetric  (or 
nearly  symmetric)  waves,  we  estimate  the  Y-wavelength  from  a  Fin* 
based  on  only  nine  gauges.  This  estimate  is  necessarily  crude. 
Consequently,  for  some  experiments  with  symmetric  waves,  we 
found  it  necessary  to  admit  two  or  three  "candidate  values"  for  Ival 
(=  Ivbl  for  a  symmetric  wave),  chosen  from  the  range  allowed  by  the 
nine-gauge  FFT. 

The  remaining  steps  of  the  algorithm  assume  that  the  wave 
pattern  is  asymmetric,  that  the  data  are  not  strictly  periodic,  and 
that  step  2  has  successfully  provided  Ival  and  Ivbl.  [For  symmetric 
waves,  P2  =  M-i.  V2  =  -  vi,  so  step  3  is  unnecessary.  Then  step  4 
simplifies,  because  d  =  b(l  -  X^).] 

Step  3:  Identify  vj  and  V2 

The  next  step  is  to  find  the  signs  of  Va  and  Vb,  and  to  determine 
which  value  goes  with  pi,  and  which  with  p2-  The  basic  idea  for  this 
step  is  the  following.  The  nine-gauge  array  lies  in  the  y-direction, 
and  the  data  from  it  determine  Ival  and  Ivbl.  If  the  array  were  rotated 
slightly  in  the  x-y  plane,  then  the  gauges  would  record  slightly  larger 
y-wavelengths  for  waves  coming  from  one  side,  and  slightly  smaller 
y-wavelengths  for  waves  coming  from  the  other.  By  comparing  the 
values  of  Ival  and  Ivbl  obtained  from  the  orginal  and  the  rotated 
gauge-arrays,  one  could  determine  whether  the  associated  waves 
came  from  the  left  or  the  right;  i.e.,  one  could  determine  the  signs  of 
Va  and  Vb- 

We  do  not  actually  move  the  gauges  to  obtain  the  second  data 
set,  but  we  achieve  approximately  the  same  effect  by  processing  data 
from  successive  gauges  at  successively  later  times.  Specifically, 
instead  of  using  the  data  from  all  nine  gauges  at  the  same  time  To, 
we  use  the  data  from  the  gauge  that  was  recorded  at  time: 

Tj  =  To  +  a  -  DAT, 

where  AT  is  the  time  between  successive  measurements.  (AT  =0.04 
seconds  for  the  2S  Hz  sampling  rate  used  in  these  experiments.) 
These  data  correspond  approximately  to  data  taken  from  a  linear 
gauge-array  whose  orientation  in  the  x-y  plane  is  rotated  from  the 
y-axis  through  an  angle 

p  =  tan-i[Vih  AT/AY],  .i7x 
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where  AY  is  the  spacing  between  gauges.  In  these  experiments, 

AY  =  2  m,  so  3  =  0.028  rad  (3=1.6®).  With  the  data  rearranged  in  this 
way,  we  repeat  step  2  to  obtain  Ivi^l  and  lvj'*'l.  The  same  procedure 
based  on  (-AT)  yields  Ivi'l  and  lvj*i.  From  all  of  these,  one  can  identify 
{Iva'I,  Iva+I)  and  {Ivb'l,  because 

Ival  -  IVa+l  -  IVa'I  -  IVal,  IVbl  -  IVb+l  -  IVb'l  -  IVbl.  (18) 

The  direction  of  effective  rotation  (3)  is  known,  so  knowing  whether 
{lva*l  <  Ival  <  Iva+I)  or  {Iva'I  >  Ival  >  Iva+I)  determines  the  sign  of  Va*,  one 
also  finds  the  sign  of  Vb  in  this  way.  [Note:  If  Ival  is  small  enough, 
then  {Va'*',  Va,  Va'}  are  not  all  of  the  same  sign,  and  no  combination  of 
measured  values  satisfies  (18).  In  this  case,  one  rinds  that  either 

IVal  +  IVa+l  «  IVa'I  -  IVal,  (X  Ival  -  IVa+l  “  IVa’I  +  IVal. 


If  both  Ival  and  Ivbl  are  small,  then  (18b)  must  be  modified  as  well.] 

The  V -values  obtained  in  this  way  also  provide  estimates  for 
(Ha.  Itb),  based  on 


Pa  ”  (Va'^  -  VaCOS  3)/sin  3* 


By  comparing  these  estimates  with  the  values  for  pi  and  P2  obtained 
in  step  1,  one  determines  (pi,vi}  and  {p2>V2}. 

Step  4:  Estimate  {b,  X,  d} 

Once  {pi,  P2.  vi,  V2)  are  known,  then  the  spatial  structure  of  the 
wave  pattern  is  fixed,  and  one  can  draw  a  period  parallelogram  in 
the  (x,  y)-plane.  Next,  one  would  like  to  choose,  among  all  KP 
solutions  of  genus  2  with  the  same  period  parallelogram,  the  one  that 
minimizes  the  rms  error. 


[f(x,y,t)  -  u(x,y,t)l^xdy 


[u(x,y,t)]2dxdy 


(19) 


where  u(x,  y,  t)  represents  the  measured  data  (normalized  as  in  (2)), 
f(x,  y,  t)  represents  a  KP  solution,  and  the  integrals  are  taken  over  a 
period  parallelogram.  In  step  6,  we  do  essentially  this. 
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Unfortunately,  apparently  has  many  local  minima  in  parameter 
space,  so  that  it  is  necessary  to  start  the  minimization  routine  close  to 
the  global  minimum.  The  purpose  of  steps  4  and  S  is  to  obtain  good 
starting  values  for  the  minimization  in  step  6. 

Let  Umax  denote  the  measured  maximum  value  of  the 
(normalized)  wave  amplitude  over  the  entire  data  set  for  an 
experiment,  and  Umin  the  measured  minimum  value.  If  the  waves 
were  exactly  periodic,  then  the  wave  amplitude  would  attain  its 
maximum  and  minimum  values  once  within  each  period 
parallelogram.  As  discussed  in  §2,  specifying  {b,  X,  d.  pi,  P2»  V2} 
specifies  a  KP  solution  of  genus  2,  up  to  a  translation.  In  particular, 
these  parameters  determine  the  maximum  (fmax)  and  minimum  (fmin) 
values  of  the  solution,  and  they  determine  vi.  Denote  by  vi.kp  the 
value  of  vi  obtained  from  a  particular  KP  solution,  and  by  vi^m  the 
value  of  Vi  measured  in  steps  2  and  3.  Step  4  is  to  minimize 

E  =  (1  -  +  (1  -  +  (1  .lmiiL)2^ 

Vl.m  “m**  nmin  (20) 

by  varying  {b,  X,  d),  subject  to  (7),  while  holding  {pi,P2.  V2}  fixed  at 
the  values  obtained  in  steps  1-3.  We  do  this  numerically,  using 
ODRPACK  as  described  by  Boggs,  Byrd,  Rogers,  and  Schnabel  (1992). 
Step  5:  Estimate  <P/  and  02 

The  KP  solution  is  now  completely  specified  except  for  a  spatial 
translation  fixed  by  the  phase  constants  and  <I>2-  These  constants 
could  be  obtained  directly  from  an  overhead  photograph  of  the 
wavefield,  such  as  those  in  Figure  1,  by  making  a  contour  plot  of  the 
KP  solution  with  the  same  horizontal  scales  as  the  photograph,  and 
then  sliding  one  over  the  other  until  the  phases  match  optimally.  In 
the  absence  of  such  a  photograph,  we  minimize  a  quantity  like  that 
in  (19),  by  varying  {0i,02}  while  holding  fixed  {b,  X,  d,  m,  p2.  V2) . 
However,  we  cannot  integrate  over  the  entire  period  parallelogram, 
as  in  (19),  but  only  over  the  nine  wave  records.  Using  ODRPACK 
again,  we  minimize 

9  /'*■ 

X  I  -  u(x,y,t)]2dx 

_ 

L  * 

[u(x,y,t)pdx  (10) 
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where  L  is  a  length  of  the  wave  record,  and  only  are  varied. 

Typically,  we  used  L  =  20  seconds. 

Step  6:  Find  final  values  of  all  eight  parameters 

The  final  step  is  to  minimize  again,  allowing  all  eight 
parameters  {b,  X,  d,  pi,  P2.  V2;  Oi,  O2}  to  vary  simultaneously.  In  this 
paper,  therefore,  the  "best"  KP  solution  for  a  given  experiment  is 
defined  to  be  the  one  that  minimizes  o^.  (In  Part  1,  we  used  a 
different  definition  for  "best"  fit  of  symmetric  waves.)  As  mentioned 
above,  has  several  local  minima  in  this  larger  space,  but  we  start 

near  what  we  believe  to  be  the  global  minimum  by  starting  at  the 
values  for  these  parameters  found  in  steps  1-5. 

Our  experience  has  been  that  this  six-step  algorithm  found  the 
best  KP  solution  for  each  of  the  fifteen  available  data  sets  for 
asymmetric  waves.  The  results  of  the  algorithm  are  discussed  in  §S. 
(It  also  converged  to  the  best  solution  for  the  data  from  sixteen  other 
experiments  on  symmetric  waves.  We  will  present  these  results 
elsewhere.) 

Our  justification  for  this  algorithm  is  that  it  works:  for  each 
experiment  on  which  we  tested  it,  the  algorithm  found  a  best  KP 
solution,  with  an  acceptable  error  (a).  However,  Boris  Dubrovin  has 
suggested  another  possible  justification  for  the  algorithm,  as  follows. 
[The  rest  of  this  paragraph  is  a  series  of  speculations,  which  we  hope 
to  verify  eventually.]  Subject  to  some  extra  conditions  like  (1),  the 
KP  equation  is  a  completely  integrable  Hamiltonian  system,  along  the 
lines  discussed  by  Dubrovin  (1991,  pp.  79-92).  Every  exact 
reduction  of  KP  to  a  finite-dimensional  system  ought  to  be 
completely  integrable  as  well.  A  KP  solution  of  genus  2  is  such  a 
reduction,  it  has  two  phases,  so  it  ought  to  correspond  to  an 
integrable  Hamiltonian  system  with  two  degrees  of  freedom. 
Therefore,  a  specification  of  the  general  KP  solution  of  genus  2  ought 
to  contain  two  "action"-type  variables,  and  two  "angles".  Steps  1-3  of 
this  algorithm  specify  the  spatial  structure,  (tii,  H2.  vi,  V2}.  The  other 
four  variables  [b,  d;  Oi,<I>2}  describe  the  dynamical  system.  Step  4, 
which  ignores  phase  information,  apparently  finds  action-type 
variables  [b,  d),  while  step  5  finds  angles  {<t>i,02}.  This  argument 
suggests  that  step  4  would  also  work  if  the  maximum  and  minimum 
values  of  the  data  set  were  replaced  by  two  conserved  integrals  of 
the  KP  equation.  In  fact,  it  might  work  better. 
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Figure  Captions 

Figure  1.  Mosaic  of  two  overhead  photographs,  showing  surface 
patterns  of  waves  in  shallow  water.  Experiment:  (a)  KP38.50204; 

(b)  KP400203-300303;  (c)  KP220203-1 10403.  Each  of  these  wave 

patterns  has  a  basic  hexagonal  template;  one  such  hexagon  is  drawn 
in  Figure  lb. 

Figure  2.  Schematic  (planform)  view  of  the  wave  basin,  showing  its 
segmented  wavemaker  (AAA/),  uniform-depth  section  (h  =  20  cm), 
and  1:30  sloping  beach  section.  The  locations  of  of  the  wave-gauges 
were  as  follows: 
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Figure  3.  Data  from  the  nine-gauge  array,  for  three  asymmetric 
wave  patterns,  including  two  shown  in  Figures  lb,  c.  Experiment: 

(a)  KP400203 -300303;  (b)  KP220203-1 10403; 

(c)  KP300102-30-803. 

Figure  4.  FFTs  (periodograms)  in  T  of  the  wave  records  shown  in 
Figure  3  for  experiment  KP220203-1 10403.  Every  FFT  shows  energy 
at  the  same  two  frequencies  (fi  =  0.36  Hz,  {2  =  0.67  Hz),  plus  their 
harmonics.  The  dominant  harmonics  have  been  identiHed  in  one  FFT. 
Figure  S.  Detailed  comparison  of  measured  wave  records  with  the 
best  KP  solution  at  the  same  locations,  for  20  seconds  of  data  from 
experiment  KP400203 -300303.  The  comparison  began  30  seconds 
into  the  experiment. 

Figure  6.  Detailed  comparison  of  measured  wave  records  with  the 
best  KP  solution,  for  20  seconds  of  data  from  experiment  KP220203- 
110403,  beginning  at  30  seconds. 

Figure  7.  For  experiment  KP220203- 110403,  contour  lines  of  the 
best  KP  solution  (identified  in  Figure  6)  are  superimposed  on  the 
photograph  of  Figure  2c.  The  experiment  and  the  KP  solution  are  the 
same  as  in  Figure  6,  but  different  data  are  compared.  The  direction 
of  propagation  of  the  entire  wave  pattern,  based  on  (8)  is  also  shown. 
Figure  8.  For  experiment  KP300102-30*803,  a  detailed  comparison 
of  measured  wave  records  with  the  best  KP  solution,  during  two 
different  time  intervals,  a)  20  seconds  of  data  beginning  at  30 
seconds;  b)  20  seconds  of  data  beginning  at  70  seconds. 

Figure  9.  Temporal  FFTs  of  the  seven  wave  records  obtained  from 
the  second  gauge  array,  for  experiment  KP220203-1 10403. 

Figure  10.  Superposition  of  the  measured  wave  records  from  the 
seven-gauge  array,  and  the  wave  pattern  predicted  by  the  KP 
solution  chosen  to  fit  the  data  from  the  nine-gauge  array,  for 
experiment  KP220203 - 1 1 0403 . 

Figure  11.  Aerial  photograph  of  waves  off  the  southern  coast  of  Long 
Island.  The  beach  is  between  Lido  Beach  and  Point  Lookout,  west  of 
Jones  Inlet.  Beyond  the  surf  zone,  the  wave  patterns  are  two- 
dimensional,  and  approximately  periodic.  They  have  flat  troughs, 
sharp  crests,  and  approximately  hexagonal  shape. 
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Figure  12.  Oblique  nonlinear  interaction  of  two  waves  in  shallow 
water,  off  the  coast  of  Oregon.  The  interaction  occurred  in  water 
about  1  m.  deep.  (Photograph  courtesy  of  T.  Toedtemeier] 

Figure  13.  Aerial  photograph  of  waves  in  shallow  water,  south  of  the 
Oregon  Inlet  on  Pea  Island,  off  the  coast  of  North  Carolina. 

[Photograph  courtesy  of  C.  Miller] 

Figure  14.  FFTs  in  Y  of  the  wave  records  for  experiment  KP220203- 
110403,  after  extension  as  described  in  the  text. 
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The  focusing  nonlinear  Schrddinger  equation  is  numerically  integrated  over  moderate  to  long  time  in¬ 
tervals.  In  certain  parameter  regimes  small  errors  on  the  order  of  roundoff  grow  rapidly  and  saturate  at 
values  comparable  to  the  main  wave.  Although  the  constants  of  motion  are  nearly  preserved,  a  serious 
phase  instability  (chaos)  develops  in  the  numerical  solutions.  The  instability  is  found  to  be  associated 
with  homoclinic  structures  and  the  underlying  mechanisms  apply  equally  well  to  many  Hamiltonian 
wave  systems. 


PACS  numbers:  03.40.Kf,  02.60.-x.  02.70. -c.  02.90.+p 

In  this  Letter  we  discuss  extensive  moderate  to  long 
time  numerical  experiments  which  we  have  carried  out  on 
the  focusing  nonlinear  Schrddinger  (NLS)  equation  with 
periodic  boundary  conditions.  The  NLS  equation  is  a 
well  known  Hamiltonian  nonlinear  wave  system  which 
arises  in  many  areas  of  physics,  and  is  special  among  such 
problems  since  a  large  class  of  solutions  can  be  computed 
via  the  inverse  scattering  transform  (1ST)  (e.g.,  [l]);  the 
NLS  equation  is  said  to  be  “integrable.”  There  are  two 
cases  of  physical  interest — the  focusing  and  defocusing 
NLS.  In  the  focusing  case,  when  periodic  boundary  con¬ 
ditions  are  imposed,  the  NLS  equation  has  complicated 
homoclinic  structures  which  under  perturbations  can  pro¬ 
duce  chaotic  dynamics  (e.g.,  (21).  The  periodic  NLS 
serves  as  a  useful  model  describing  unstable  wave  phe¬ 
nomena  (e.g.,  instability  in  deep  water  waves)  and  has 
been  the  subject  of  numerical  simulations  and  laboratory 
experiments  (e.g.,  (3]). 

In  our  investigations,  we  employ  two  numerical 
schemes  which  have  been  used  extensively  and  effectively 
by  researchers  studying  the  NLS  equation:  (a)  the  in¬ 
tegrable  discrete  NLS  (IDNLS)  equation  (e.g.,  (2,41) 
and  (b)  the  Fourier  split-step  (FSS)  algorithm  (e.g.,  [5l). 
The  IDNLS  equation  is  an  integrable  differential-dif¬ 
ference  equation  [4]  and  is  implemented  using  a  high  or¬ 
der  time  discretization.  The  FSS  algorithm,  although  not 
integrable,  preserves  the  underlying  symplectic  structure 
of  the  NLS  equation  and,  as  such,  is  in  the  class  of  sym¬ 
plectic  integrators  which  have  been  used  as  a  means  of 
tracking  the  long  time  behavior  of  Hamiltonian  systems 


(e.g.,  [6]  and  references  therein).  We  use  two  numerical 
schemes  to  demonstrate  that  the  results  obtained  are  due 
to  the  extreme  sensitivity  of  the  periodic  focusing  NLS 
equation  (in  the  parameter  regime  described  below)  and 
not  the  particular  details  of  the  numerical  schemes  em¬ 
ployed.  To  be  brief  we  mainly  discuss  the  calculations  of 
IDNLS.  The  FSS  algorithm  yields  analogous  results. 

In  earlier  work  we  have  shown  that  initial  data  which 
are  nearby  low  dimensional  "homoclinic  manifolds” 
trigger  numerically  induced  joint  spatial  and  temporal 
chaos  in  nonintegrable  numerical  schemes  at  intermediate 
values  of  the  mesh  size  [2].  This  chaos  disappears  as  the 
mesh  is  refined.  In  this  Letter  we  concentrate  on  a  more 
troubling  aspect  of  numerically  induced  chaos.  We  show 
that  temporal  instabilities  and  chaos  can  be  easily  excited 
by  very  small  perturbations — on  the  order  of  roundoff. 
Although  our  discussion  centers  on  the  NLS  equation,  we 
have  observed  analogous  behavior  in  other  problems  such 
as  the  sine-Gordon  and  modified  Korteweg-de  Vries 
(KdV)  equations.  We  believe  that  similar  results  will  be 
found  in  many  other  Hamiltonian  systems.  The  NLS 
equation  is  an  excellent  paradigm  system  to  study  since 
we  have  a  great  deal  of  analytical  knowledge  about  this 
equation,  and  it  is  reasonably  straightforward  to  com¬ 
pute. 

We  begin  by  summarizing  our  main  observations. 

(I)  Tiny  numerical  errors  (i.e.,  10“'*)  grow  rapidly, 
eventually  saturate,  but  significantly  alter  the  solution 
after  moderate  times.  For  example,  spatially  even  initial 
values  must  evolve  in  an  even  manner.  However,  without 
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imposing  evenness  as  a  separate  constraint,  we  find  that 
an  odd  component  is  excited  and  develops  into  a  size 
comparable  with  the  “main  wave."  Alternatively,  calcu¬ 
lating  in  mathematically  but  not  computationally  equiv¬ 
alent  ways  (e.g.,  simply  by  reordering  terms  in  the  equa¬ 
tion)  shows  that  small  errors  on  the  order  of  roundolT 
grow  rapidly  and  eventually  destroy  the  “true"  solution. 

(2)  In  the  Hamiltonian  framework,  the  periodic  spec¬ 
trum  (see  below)  provides  the  “actions"  in  an  action- 
angle  description  and  they  are  constants  of  the  motion. 
We  verify  that  these  constants  are  preserved  by  the  nu¬ 
merical  scheme  to  very  high  precision. 

(3)  The  sensitivity  is  due  to  the  proximity  to  the  under¬ 
lying  homoclinic  orbits  of  the  NLS  equation.  It  (urns  out 
that  for  the  range  of  initial  values  we  have  chosen  the 
splitting  distance  between  suitable  complex  eigenvalues  of 
the  periodic  spectrum  is  extremely  small.  Such  initial 
data  occur  naturally  and  can  be  associated  with  the  evo¬ 
lution  of  nearly  elementary  plane  waves  as  well  as  waves 
which  evolve  into  complicated  nonlinear  states. 

(4)  After  the  growth  of  the  small  errors,  the  corre¬ 
sponding  state  evolves,  but  not  in  a  quasiperiodic  tem¬ 
poral  manner  as  might  have  been  expected  based  on  the 
integrability  of  NLS.  The  power  spectrum  of  the  associ¬ 
ated  evolution  has  a  slowly  decaying  “tail"  whose  size 
grows  as  the  number  of  linearly  unstable  modes  <i.e.,  the 
suitable  complex  eigenvalues)  increases. 

(5)  Analogous  solutions  of  the  defocusing  NLS  are 
stable.  They  are  net  sensitive,  and  we  have  no  difficulty 
in  computing  the  solution  over  long  times,  using  either 
IDNLS  or  FSS.  The  power  spectrum  remains  compact 
and,  computationally,  the  solution  is  quasiperiodic  in 
time.  The  inverse  scattering  theory  of  the  focusing  NLS 
equation  with  periodic  boundary  values  is  significantly 
different  from  the  defocusing  NLS  or  the  KdV  equation 
[7].  Defocusing  NLS  and  KdV  are  analogous  to  finite  di¬ 
mensional  Hamiltonian  systems  whose  underlying  geo¬ 
metric  structures  are  compact  tori.  In  these  cases, 
without  homoclinic  structures,  symplectic  integrators  are 
effective. 

We  consider  the  NLS  equation  in  the  form 

iu,  +  uxx  +  2s\u\^u“‘0 ,  (I) 

where  s“±I  (-1-1  focusing,  —I  defocusing),  with 
periodic  boundary  conditions  «(jf  —  L/2)  “u(x -F/./2) 
and  initial  values  u(x,0)  "gfx).  In  the  focusing  case, 
the  plane  wave  solution,  mo(x,/) “flexp(2i |a|^r),  is 
linearly  unstable.  When  considering  small  perturbations 
of  the  form  «(x,r) “«o(jc,/)Il +f(x,/)l,  one  can  estab¬ 
lish  that  the  perturbation  is  separated  into  solutions  of 
the  form  f(x,i)a:exp(ifi„x)exp(<T„i),  where 
and  the  growth  rate  a„  is  given  by  a„  —/inV^loP— /“«• 
Consequently,  the  solution  is  unstable  provided 
<  4|a|  .  The  number  of  unstable  modes  is  the  largest  in¬ 
teger  M  satisfying  0<  M  <\a\Lln.  The  fastest  growth 
rate  corresponds  to  a  value  n(ji„)  for  which  a„  is  maxi¬ 


mal  (i.e.,  closest  to  2|a|^). 

To  solve  the  NLS  via  1ST  we  use  the  associated  linear 
scattering  problem 

J’lx  — /X,i'i-Fm(x,/)i’2  ,  (2) 

l'2x“i^l'2~SU*(x,l)l'\  ,  (3) 

where  u  *  is  the  complex  conjugate  of  u  and  the  parame¬ 
ter  X  is  an  eigenvalue  once  boundary  values  are  specified. 
The  spectrum  of  this  linear  operator  is  invariant  under 
the  NLS  flow  and  the  periodic  or  antiperiodic  eigenvalues 
provide  sufficient  invariants  to  establish  the  in¬ 
tegrability  of  the  N  LS.  The  conserved  integrals  S-l/2Ti 
xrfx“*C/,i“l,2,.. .  such  as  the  norm  Ti—luP,  the 
momentum  and  the  Hamiltonian  T's'-luxP 

-l-|u|^  are  related  to  the  above  eigenvalues. 

The  IDNLS  equation  is  given  by  [2,4] 

/M„-»-(u,  +  i+«,-i  -2u«)//i^-l-s|u„|^(u„  +  i-t-u„-i)  -0 

(4) 

with  u„(r)  ■■u(nfi,r),  Uii-w/2(/)“«(i+w/2(t),  n"0, 

I, ...  .A'  {N  even),  h  “L/N,  s  “  ±  1.  Since  IDNLS  is 
integrable,  it  possesses  the  special  properties  inherent  in 
such  systems.  As  A— *0  (4)  converges  to  (1)  with  accu¬ 
racy  O(A^).  In  our  numerical  simulations  the  time  in¬ 
tegration  of  this  system  of  ordinary  differential  equations 
is  performed  by  the  adaptive  Runge-Kutta-Merson  rou¬ 
tine  (d02BBF)  in  the  NAG  (numerical  algorithms  group) 
software  library.  We  specify  an  extremely  high  tolerance 
to  ensure  that  the  results  are  not  dependent  upon  the  time 
integration.  For  a  description  of  the  FSS  algorithm  see, 
for  example,  (5].  We  only  remark  that  in  each  partial 
step  of  the  FSS  algorithm,  the  canonical  symplectic 
structure  of  the  NLS  equation  is  preserved,  hence  the 
FSS  algorithm  is  a  symplectic  integrator. 

We  first  discuss  our  calculations  associated  with  the 
initial  data  g(x,e)’“fl(I-l-fcos6px),  where  a"*1.5,  p 
““2k/L,  L  and  f  ■•10“'^  From  the  linear  theory 

discussed  above,  we  see  that  g(x,f )  contains  eight  unsta¬ 
ble  modes.  In  Fig.  1  we  plot  the  results  of  IDNLS  where 
we  use  512  grid  points.  In  Fig.  1(a)  we  plot  the  growth 
of  asymmetry  [from  (1)  note  that  g(x,e)  is  an  even  func¬ 
tion  and  should  have  evolved  in  an  even  manner].  The 
asymmetry  is  measured  by  logigl/*— /-4I  where  ft  is  the 
sixth  Fourier  harmonic  [u  “J^mCxpO'pmx)].  The  asym¬ 
metry  grows  exponentially  fast;  for  a  while  the  rate  of 
growth  is  consistent  with  that  predicted  by  the  maximum 
growth  rate  of  the  linearly  unstable  plane  wave  (i.c., 
—4.5).  Note  that  the  asymmetry  reaches  0(1)  be¬ 
fore  it  saturates.  In  Fig.  1(b)  an  averaged  amplitude  of 
the  Fourier  transform  is  depicted.  Let  uj(f:)  be  the 
Fourier  transform  of  uj(l)  luj(l)  is  evaluated  at  every 
A/ ’■0.2];  we  plot  the  average:  /,,e(k) 
vs  k.  We  note  that  the  Fc  ricr  spectrum  has  a  slowly  de¬ 
caying  asymptotic  “tail"  which  is  inconsistent  with  the 
expectation  of  a  quasiperiodic  temporal  flow.  Note  espe- 
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tim« 

FIG.  I.  (a)  The  asymmetry  for  IDNLS  with  "10”'*). 
(b)  The  power  spectrum  for  IDNLS  with  “10”'*).  (c) 
The  growth  in  the  difference  of  solutions  for  IDNLS  with 
gfx.f-IO"'*). 

cially  the  comparison  with  the  Fourier  spectrum  of  the 
defocusing  NLS  equation  (Fig.  4).  In  Fig.  I  (c)  we  plot 
the  logarithm  of  the  averaged  difference  of  two  solutions; 
logioZj^olMy(f)  —“y(f)|/((V  +  I )  where  ujit)  is  the  solu¬ 
tion  calculated  on  i  —  Ll2,Ll2)  and  ujit)  is  the  solution 
calculated  on  (0,L/2)  with  symmetry  imposed,  i.e., 
u  -„(/)  “u„(l)  (i.e.,  we  calculated  on  half  the  lattice).  In 
these  calculations  the  “dominant”  complex  eigenvalues  in 
the  upper  half  plane  were  nearly  isospectral.  The  relative 
change  of  the  first  five  conserved  quantities  {C,}/— |  was 


FIG.  2.  The  growth  in  the  difference  of  solutions  for  IDNLS 
with  g(jt,e“IO”'). 


less  than  0.01%.  As  we  vary  e  in  g(x,f)  we  find  similar 
results;  although  for  significantly  larger  f  (e.g.,  e"IO~') 
it  takes  longer  before  the  asymmetric  perturbations  “or¬ 
ganize”  into  rapid  growth. 

As  a  crude  model  of  roundoff  effects  we  add  to  the 
solution  [e.g.,  g(x,e)  with  e~l0~']  a  small  random  per¬ 
turbation  at  each  time  step:  ftr(x,l),  where  0  <  r  <  I  (r 
is  obtained  via  a  random  number  generator)  for  a  range 
of  fi:  10"*,  10"’ . 10"'*.  The  main  features  ob¬ 

served  before  hold:  rapid  growth  of  asymmetry,  a  slowly 
decaying  tail  in  the  power  spectrum,  etc. 

Even  if  we  impose  evenness  as  a  constraint,  errors  on 
the  order  of  roundoff  grow  rapidly.  For  example,  consid¬ 
er  IDNLS  with  the  initial  values  g(x,f)  with  e"10~' 
and  u-„(t)^u„it)  imposed.  We  first  calculate  IDNLS 
in  precisely  the  way  Eq.  (4)  is  specified  for  r— 0-100 
— call  this  solution  u^'^  Next  we  calculate  (4)  only  by 
distributing  the  nonlinear  terms,  i.e.,  using  5|u„|*Un.|.i 
Call  this  solution  The  logarithm  of 
the  averaged  difference  between  the  solutions, 
Zi'^SlogiolKy'Hf)  —  «/^Hf)|/(Af/2-H ),  grows  rapidly 
(Fig.  2)  and,  as  was  the  case  with  asymmetry,  the  solu¬ 
tions  become  drastically  different.  We  note  that  the 
power  spectra  of  each  is  broad  banded  and  different  from 
one  another  as  well  (not  plotted  here). 

Strongly  nonlinear  states  are  unstable  and  experience 
the  same  loss  of  predictability  associated  with  the  above 
phenomena.  For  initial  data  gix,e)  witii  f“*IO"'  we 
calculate  the  solution  via  IDNLS  for  t  “0-25.  At  /  “25 
we  are  into  the  fully  nonlinear  regime  and  using  these 
values  as  initial  data,  we  then  calculate  the  two  solutions 
(above):  from  /“25  to  /“225.  Again 

(Fig.  3)  the  solutions  deviate  significantly,  although  the 
growth  rate  of  the  difference  is  smaller  than  in  the  previ¬ 
ous  case.  Furthermore,  imposing  a  small  random  pertur¬ 
bation  at  each  time  step  (for  f  >  25),  we  find  a  more  rap¬ 
id  growth  of  the  deviation  to  an  0(1 )  value. 

For  defocusing  NLS  all  the  difficulties  encountered  in 
the  focusing  case  disappear.  The  power  spectrum  of  a 
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FIG.  3.  The  growth  in  the  dilTerence  of  solutions  for  IDNLS  FIG.  4.  The  power  spectrum  for  dcfocusing  IDNLS  with 
for  the  nonlinear  regime.  g(*,e "10“'). 


typical  solution  [compare  with  Fig.  1(b)]  is  plotted  in 
Fig.  4  where  we  take  g(x,^),e■■IO”'.  We  observe  no 
growth  of  asymmetry  and  both  algorithms  IDNLS  and 
FSS  show  no  growth  in  the  difference  between  solutions 
and  Recall  that  the  periodic  dcfocusing  NLS 
does  not  have  homoclinic  structures.  Small  errors  do  not 
grow  rapidly  and  IDNLS  and  FSS  are  accurate  predic¬ 
tors  of  the  long  time  behavior  of  the  solution. 

Since  we  are  studying  initial  data  which  are  small  devi¬ 
ations  from  a  constant,  we  are  able  to  calculate  the 
periodic  or  antiperiodic  eigenvalues  of  2-3  via  perturba¬ 
tion  analysis  (see  e.g.,  [8]).  For  the  initial  data  g(x,e) 
■■fl(l  +  fcos6/ix)  with  (  small,  at  leading  order  there  are 
eight  complex  eigenvalues  which  are  double  points  (i.e., 
eigenvalues  of  multiplicity  2).  Each  of  these  double 
points  corresponds  to  an  unstable  mode.  We  find  that  un¬ 
der  perturbation  the  double  point  corresponding  to  the 
perturbation  6fi  is  split  by  Oie).  The  other  complex  dou¬ 
ble  points,  while  translating  by  Oif^)  nevertheless  have  a 
splitting  distance  which  is  smaller  than  any  power  of  (. 
The  splitting  distance  of  the  remaining  complex  eigenval¬ 
ues  is  beyond  all  orders  in  f.  Consequently,  we  are  ’’ex¬ 
ponentially  close”  to  the  homoclinic  manifolds,  and  small 
deviations  on  the  order  of  roundoff  can  lead  to  homoclinic 
crossings  and  the  observed  chaos.  We  have  extended  our 
perturbation  analysis  to  predict  how  modifications  in  the 
initial  values  of  earlier  numerical  simulations  [3]  could 
lead  to  homoclinic  crossings  under  tiny  perturbations 
such  as  roundoff.  We  believe  that  this  mechanism  can 
also  explain  similar  difficulties  observed,  but  not  under¬ 
stood,  in  other  Hamiltonian  equations  [9]. 

The  extreme  proximity  of  these  data  values  to  the  un¬ 
perturbed  homoclinic  manifolds  prevents  us  from  numeri¬ 
cally  calculating  the  “true”  solution  after  moderate  time. 
The  calculations  are  influenced  by  minuscule  errors;  com¬ 
putationally  the  NLS  equation  is  “effectively  chaotit;”  in 
this  range  of  parameter  space  despite  it  being  “inte- 
grable.”  Since  NLS  is  so  special,  we  expect  an  underly¬ 
ing  structure  to  the  temporal  disorder  which  may  be 
determined  by  appropriate  statistical  studies.  Numerical 
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simulations  which  test  the  stability  of  small  odd  perturba¬ 
tions,  rearrangements  of  terms,  complex  nonlinear  states, 
such  as  those  mentioned  in  this  Letter,  can  demonstrate 
the  existence  of  serious  underlying  numerical  instabilities 
even  when  known  constants  of  the  motion  are  seemingly 
well  preserved.  A  more  complete  account  of  the  work  de¬ 
scribed  here  will  be  published  in  a  future  paper. 
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